EDUC 152. Intro to quantitative research in education: Regression analysis
Statistical inference
1 Introduction
Statistical inference
The theory, methods, and practice of forming judgments about the parameters of a population and the reliability of statistical relationships, typically on the basis of random sampling
credit: Google
1.1 Libraries, data, functions
Purpose of this section:
- load all libraries/packages we will use
- loads and creates IPEDS data frame (population)
- creates data frame of generated variables (population)
- creates sample versions of the IPEDS and gnerated data frames
- create all functions we will use
- credit to data scientist extraordinaire Crystal Han for creating most of these functions!
Why do all of this in one place?
- so you don’t have to run code chunks scattered randomly throughout the lecture
- dont worry about understanding all code, especially for functions we create
# Libraries
#install.packages('tidyverse') # if you haven't installed already
#install.packages('labelled') # if you haven't installed already
#install.packages('patchwork') # if you haven't installed already
library(tidyverse) # load tidyverse package
library(labelled) # load labelled package package
library(patchwork)
##########
########## IPEDS
##########
# Load ipeds dataset from course website url
load(file = url('https://github.com/anyone-can-cook/educ152/raw/main/data/ipeds/output_data/panel_data.RData'))
# Create ipeds data frame with fewer variables/observations
df_ipeds_pop <- panel_data %>%
# keep data from fall 2019
filter(year == 2019) %>%
# which universities to keep:
# 2015 carnegie classification: keep research universities (15,16,17) and master's universities (18,19,20)
filter(c15basic %in% c(15,16,17,18,19,20)) %>%
# which variables to keep
select(instnm,unitid,opeid6,opeid,control,c15basic,stabbr,city,zip,locale,obereg, # basic institutional characteristics
tuition6,fee6,tuition7,fee7, # avg tuition and fees for full-time grad, in-state and out-of-state
isprof3,ispfee3,osprof3,ospfee3, # avg tuition and fees for MD, in-state and out-of-state
isprof9,ispfee9,osprof9,ospfee9, # avg tuition and fees for Law, in-state and out-of-state
chg4ay3,chg7ay3,chg8ay3) %>% # [undergraduate] books+supplies; off-campus (not with family) room and board; off-campus (not with family) other expenses
# rename variables; syntax <new_name> = <old_name>
rename(region = obereg, # revion
tuit_grad_res = tuition6, fee_grad_res = fee6, tuit_grad_nres = tuition7, fee_grad_nres = fee7, # grad
tuit_md_res = isprof3, fee_md_res = ispfee3, tuit_md_nres = osprof3, fee_md_nres = ospfee3, # md
tuit_law_res = isprof9, fee_law_res = ispfee9, tuit_law_nres = osprof9, fee_law_nres = ospfee9, # law
books_supplies = chg4ay3, roomboard_off = chg7ay3, oth_expense_off = chg8ay3) %>% # [undergraduate] expenses
# create measures of tuition+fees
mutate(
tuitfee_grad_res = tuit_grad_res + fee_grad_res, # graduate, state resident
tuitfee_grad_nres = tuit_grad_nres + fee_grad_nres, # graduate, non-resident
tuitfee_md_res = tuit_md_res + fee_md_res, # MD, state resident
tuitfee_md_nres = tuit_md_nres + fee_md_nres, # MD, non-resident
tuitfee_law_res = tuit_law_res + fee_law_res, # Law, state resident
tuitfee_law_nres = tuit_law_nres + fee_law_nres) %>% # Law, non-resident
# create measures of cost-of-attendance (COA) as the sum of tuition, fees, book, living expenses
mutate(
coa_grad_res = tuit_grad_res + fee_grad_res + books_supplies + roomboard_off + oth_expense_off, # graduate, state resident
coa_grad_nres = tuit_grad_nres + fee_grad_nres + books_supplies + roomboard_off + oth_expense_off, # graduate, non-resident
coa_md_res = tuit_md_res + fee_md_res + books_supplies + roomboard_off + oth_expense_off, # MD, state resident
coa_md_nres = tuit_md_nres + fee_md_nres + books_supplies + roomboard_off + oth_expense_off, # MD, non-resident
coa_law_res = tuit_law_res + fee_law_res + books_supplies + roomboard_off + oth_expense_off, # Law, state resident
coa_law_nres = tuit_law_nres + fee_law_nres + books_supplies + roomboard_off + oth_expense_off) %>% # Law, non-resident
# keep only observations that have non-missing values for the variable coa_grad_res
# this does cause us to lose some interesting universities, but doing this will eliminate some needless complications with respect to learning core concepts about statistical inference
filter(!is.na(coa_grad_res))
# Add variable labels to the tuit+fees variables and coa variables
# tuition + fees variables
var_label(df_ipeds_pop[['tuitfee_grad_res']]) <- 'graduate, full-time, resident; avg tuition + required fees'
var_label(df_ipeds_pop[['tuitfee_grad_nres']]) <- 'graduate, full-time, non-resident; avg tuition + required fees'
var_label(df_ipeds_pop[['tuitfee_md_res']]) <- 'MD, full-time, state resident; avg tuition + required fees'
var_label(df_ipeds_pop[['tuitfee_md_nres']]) <- 'MD, full-time, non-resident; avg tuition + required fees'
var_label(df_ipeds_pop[['tuitfee_law_res']]) <- 'Law, full-time, state resident; avg tuition + required fees'
var_label(df_ipeds_pop[['tuitfee_law_nres']]) <- 'Law, full-time, non-resident; avg tuition + required fees'
# COA variables
var_label(df_ipeds_pop[['coa_grad_res']]) <- 'graduate, full-time, state resident COA; == tuition + fees + (ug) books/supplies + (ug) off-campus room and board + (ug) off-campus other expenses'
var_label(df_ipeds_pop[['coa_grad_nres']]) <- 'graduate, full-time, non-resident COA; == tuition + fees + (ug) books/supplies + (ug) off-campus room and board + (ug) off-campus other expenses'
var_label(df_ipeds_pop[['coa_md_res']]) <- 'MD, full-time, state resident COA; == tuition + fees + (ug) books/supplies + (ug) off-campus room and board + (ug) off-campus other expenses'
var_label(df_ipeds_pop[['coa_md_nres']]) <- 'MD, full-time, non-resident COA; == tuition + fees + (ug) books/supplies + (ug) off-campus room and board + (ug) off-campus other expenses'
var_label(df_ipeds_pop[['coa_law_res']]) <- 'Law, full-time, state resident COA; == tuition + fees + (ug) books/supplies + (ug) off-campus room and board + (ug) off-campus other expenses'
var_label(df_ipeds_pop[['coa_law_nres']]) <- 'Law, full-time, non-resident COA; == tuition + fees + (ug) books/supplies + (ug) off-campus room and board + (ug) off-campus other expenses'
df_ipeds_pop %>% head()#> # A tibble: 6 x 38
#> instnm unitid opeid6 opeid control c15basic stabbr city zip locale region tuit_grad_res fee_grad_res tuit_grad_nres fee_grad_nres tuit_md_res fee_md_res tuit_md_nres fee_md_nres tuit_law_res fee_law_res tuit_law_nres fee_law_nres books_supplies roomboard_off oth_expense_off tuitfee_grad_res tuitfee_grad_nres tuitfee_md_res tuitfee_md_nres tuitfee_law_res tuitfee_law_nres coa_grad_res coa_grad_nres coa_md_res coa_md_nres coa_law_res coa_law_nres
#> <chr> <dbl> <chr> <chr> <dbl+lbl> <dbl+lbl> <chr+lbl> <chr> <chr> <dbl+lbl> <dbl+lbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
#> 1 Alabama A & M University 100654 001002 00100200 1 [Public] 18 [Master^s Colleges & Universities: Larger Programs] AL [Alabama] Normal 35762 12 [City: Midsize] 5 [Southeast AL AR FL GA KY LA MS NC SC TN VA WV] 10128 1414 20160 1414 NA NA NA NA NA NA NA NA 1600 9240 3090 11542 21574 NA NA NA NA 25472 35504 NA NA NA NA
#> 2 University of Alabama at Birmingham 100663 001052 00105200 1 [Public] 15 [Doctoral Universities: Highest Research Activity] AL [Alabama] Birmingham 35294-0110 12 [City: Midsize] 5 [Southeast AL AR FL GA KY LA MS NC SC TN VA WV] 8100 0 19188 0 28978 0 62714 0 NA NA NA NA 1200 12307 5555 8100 19188 28978 62714 NA NA 27162 38250 48040 81776 NA NA
#> 3 Amridge University 100690 025034 02503400 2 [Private not-for-profit] 20 [Master^s Colleges & Universities: Small Programs] AL [Alabama] Montgomery 36117-3553 12 [City: Midsize] 5 [Southeast AL AR FL GA KY LA MS NC SC TN VA WV] 11700 1300 11700 1300 NA NA NA NA NA NA NA NA 900 9600 1600 13000 13000 NA NA NA NA 25100 25100 NA NA NA NA
#> 4 University of Alabama in Huntsville 100706 001055 00105500 1 [Public] 16 [Doctoral Universities: Higher Research Activity] AL [Alabama] Huntsville 35899 12 [City: Midsize] 5 [Southeast AL AR FL GA KY LA MS NC SC TN VA WV] 10632 826 24430 826 NA NA NA NA NA NA NA NA 2120 10400 3994 11458 25256 NA NA NA NA 27972 41770 NA NA NA NA
#> 5 Alabama State University 100724 001005 00100500 1 [Public] 19 [Master^s Colleges & Universities: Medium Programs] AL [Alabama] Montgomery 36104-0271 12 [City: Midsize] 5 [Southeast AL AR FL GA KY LA MS NC SC TN VA WV] 7416 2740 14832 2740 NA NA NA NA NA NA NA NA 1600 7320 4228 10156 17572 NA NA NA NA 23304 30720 NA NA NA NA
#> 6 The University of Alabama 100751 001051 00105100 1 [Public] 16 [Doctoral Universities: Higher Research Activity] AL [Alabama] Tuscaloosa 35487-0100 12 [City: Midsize] 5 [Southeast AL AR FL GA KY LA MS NC SC TN VA WV] 10780 0 30250 0 28978 0 62714 0 23610 0 43060 0 1000 13636 4600 10780 30250 28978 62714 23610 43060 30016 49486 48214 81950 42846 62296##########
########## Create data frame of generated variables, with each variable meant to represent the entire population
##########
num_obs <- 10000
# Generate normal distribution w/ custom mean and sd
set.seed(124)
norm_dist <- rnorm(n = num_obs, mean = 50, sd = 5)
# Generate right-skewed distribution
set.seed(124)
rskew_dist <- rbeta(n = num_obs, shape1 = 2, shape2 = 5)
# Generate left-skewed distribution
set.seed(124)
lskew_dist <- rbeta(n = num_obs, shape1 = 5, shape2 = 2)
# Generate standard normal distribution (default is mean = 0 and sd = 1)
set.seed(124)
stdnorm_dist <- rnorm(n = num_obs, mean = 0, sd = 1) # equivalent to rnorm(10)
# Create dataframe
df_generated_pop <- data.frame(norm_dist, rskew_dist, lskew_dist, stdnorm_dist) %>% as.tibble()
# drop individual objects associated with each variable
rm(norm_dist,rskew_dist,lskew_dist,stdnorm_dist)
rm(num_obs)
##########
########## Create sample versions of generated population data frame and IPEDS population data frame
##########
# create sample version of our generated data
set.seed(124) # set seed so that everyone ends up with the same random sample
df_generated_sample <- df_generated_pop %>% sample_n(size = 200)
df_generated_sample %>% head()#> # A tibble: 6 x 4
#> norm_dist rskew_dist lskew_dist stdnorm_dist
#> <dbl> <dbl> <dbl> <dbl>
#> 1 42.7 0.345 0.655 -1.46
#> 2 50.2 0.314 0.686 0.0488
#> 3 61.3 0.0938 0.906 2.26
#> 4 45.5 0.0558 0.944 -0.905
#> 5 44.7 0.0745 0.926 -1.05
#> 6 48.0 0.301 0.699 -0.402# create sample version of our ipeds data
set.seed(124) # set seed so that everyone ends up with the same random sample
df_ipeds_sample <- df_ipeds_pop %>% sample_n(size = 200)
# compare mean of coa_grad_res between population and sample
mean(df_ipeds_pop$coa_grad_res, na.rm = TRUE)#> [1] 30318.74 mean(df_ipeds_sample$coa_grad_res, na.rm = TRUE)#> [1] 30002.74##########
# Write function to get the sampling distribution from a variable (defaults equal 1000 samples of size 200)
##########
get_sampling_distribution <- function(data_vec, num_samples = 1000, sample_size = 200) {
sample_means <- vector(mode = 'numeric', num_samples)
for (i in 1:length(sample_means)) {
samp <- sample(data_vec, sample_size)
sample_means[[i]] <- mean(samp, na.rm = T)
}
sample_means
}
##########
# Create function to generate plots of variable distributions
##########
plot_distribution <- function(data_df, data_var, group_var = NULL, group_cat = NULL, show_group_hist = F, sampling_dist = F, plot_title = '') {
# Prep dataframe
two_pop <- !is.null(group_var)
two_dist <- !is.null(group_var) && !sampling_dist
data_df[[data_var]] <- unclass(data_df[[data_var]]) # unclass haven_labelled
if (two_pop) {
data_df[[group_var]] <- unclass(data_df[[group_var]])
}
data_df <- data_df %>% filter(!is.na(get(data_var))) # remove NA rows
# If group_cat not provided, use 2 values from group_var
if (two_pop && is.null(group_cat)) {
group_cat <- sort(unique(na.omit(data_df[[group_var]])))[1:2]
}
# Create population vector(s)
if (!two_pop) { # single population
data_vec1 <- data_df[[data_var]]
if (sampling_dist) { # sampling distribution
data_vec1 <- get_sampling_distribution(data_vec1)
}
} else { # two populations
data_vec1 <- (data_df %>% filter(get(group_var) == group_cat[[1]]))[[data_var]]
data_vec2 <- (data_df %>% filter(get(group_var) == group_cat[[2]]))[[data_var]]
if (sampling_dist) { # sampling distribution
data_vec1_samp <- get_sampling_distribution(data_vec1)
data_vec2_samp <- get_sampling_distribution(data_vec2)
data_vec1 <- data_vec1_samp - data_vec2_samp
}
}
# Create statistics dataframe
if (!two_dist) {
lines_vec <- c('dotted')
stats_vec <- c(mean(data_vec1), median(data_vec1))
legend_title <- 'Statistics'
if (two_pop && sampling_dist) {
legend_title <- paste0('Statistics\n(', group_var, '=', group_cat[[1]], ' - ', group_var, '=', group_cat[[2]], ')')
}
} else {
lines_vec <- c('dotted', 'dotdash')
stats_vec <- c(mean(data_vec1), median(data_vec1), mean(data_vec2), median(data_vec2))
legend_title <- paste0('Statistics\n(', group_var, '=', group_cat[[1]], ' vs. ', group_var, '=', group_cat[[2]], ')')
}
stats_df <- data.frame(
pop = rep(lines_vec, each = 2),
stat = rep(c('blue', 'red'), times = if_else(two_dist, 2, 1)),
val = stats_vec
)
stats_df$pop <- factor(stats_df$pop, levels = c('dotted', 'dotdash'))
# Legend text
legend_text <- c(paste('Mean:', round(mean(data_vec1), 2),
'\nStd Dev:', round(sd(data_vec1), 2)),
paste('Median:', round(median(data_vec1), 2)))
if (two_dist) {
legend_text <- c(legend_text,
paste('Mean:', round(mean(data_vec2), 2),
'\nStd Dev:', round(sd(data_vec2), 2)),
paste('Median:', round(median(data_vec2), 2)))
}
# Plot distribution(s)
p <- ggplot() +
ggtitle(plot_title) + xlab('') + ylab('') +
geom_density(aes(x = data_vec1), alpha = 0.8)
if (!two_dist || show_group_hist) { # show histogram only if 1 pop or show_group_hist is TRUE
p <- p +
geom_histogram(aes(x = data_vec1, y = ..density..), alpha = 0.4, position = 'identity')
}
if (two_dist) {
p <- p +
geom_density(aes(x = data_vec2), alpha = 0.8)
if (show_group_hist) { # show histogram only if show_group_hist is TRUE
p <- p +
geom_histogram(aes(x = data_vec2, y = ..density..), alpha = 0.4, position = 'identity', fill = 'wheat4')
}
}
p <- p +
geom_vline(data = stats_df,
aes(xintercept = val, color = interaction(stat, pop), linetype = interaction(stat, pop)),
size = 0.6, alpha = 0.8) +
scale_color_manual(name = legend_title,
labels = legend_text,
values = as.character(stats_df$stat)) +
scale_linetype_manual(name = legend_title,
labels = legend_text,
values = as.character(stats_df$pop)) +
theme(plot.title = element_text(size = 10, face = 'bold', hjust = 0.5),
legend.title = element_text(size = 9, face = 'bold'),
legend.text = element_text(size = 8)) +
guides(col = guide_legend(ncol = if_else(two_dist, 2, 1)))
p
}
##########
# Write Function to generate sampling distribution (with t-test value) assuming null hypothesis is correct
##########
# Function to generate t-distribution plot
plot_t_distribution <- function(data_df, data_var, group_var = NULL, group_cat = NULL, mu = 0, alpha = 0.05, alternative = 'two.sided', plot_title = '', shade_rejection = T, shade_pval = T, stacked = F) {
# Prep dataframe
two_pop <- !is.null(group_var)
data_df[[data_var]] <- unclass(data_df[[data_var]]) # unclass haven_labelled
if (two_pop) {
data_df[[group_var]] <- unclass(data_df[[group_var]])
}
data_df <- data_df %>% filter(!is.na(get(data_var))) # remove NA rows
# If group_cat not provided, use 2 values from group_var
if (two_pop && is.null(group_cat)) {
group_cat <- sort(unique(na.omit(data_df[[group_var]])))[1:2]
}
# Calculate stats
if (!two_pop) { # single sample
data_vec <- data_df[[data_var]]
# Calculate t-statistics
sample_size <- length(data_vec)
deg_freedom <- sample_size - 1
xbar <- mean(data_vec)
s <- sd(data_vec)
std_err <- s / sqrt(sample_size)
t <- (xbar - mu) / std_err
} else { # two samples
data_vec1 <- (data_df %>% filter(get(group_var) == group_cat[[1]]))[[data_var]]
data_vec2 <- (data_df %>% filter(get(group_var) == group_cat[[2]]))[[data_var]]
# Calculate t-statistics
xbar1 <- mean(data_vec1)
xbar2 <- mean(data_vec2)
s1 <- sd(data_vec1)
s2 <- sd(data_vec2)
n1 <- length(data_vec1)
n2 <- length(data_vec2)
deg_freedom <- (s1**2/n1 + s2**2/n2)**2 / ((s1**2/n1)**2/(n1-1) + (s2**2/n2)**2/(n2-1))
std_err <- sqrt(s1**2/n1 + s2**2/n2)
t <- (xbar1 - xbar2) / std_err
}
# Calculate critical value and p-value
if (alternative == 'less') { # left-tailed
cv_lower <- qt(p = alpha, df = deg_freedom, lower.tail = T)
cv_legend <- round(cv_lower, 2)
cv_legend2 <- round(cv_lower * std_err + mu, 2)
pval <- round(pt(q = t, df = deg_freedom, lower.tail = T), 4)
} else if (alternative == 'greater') { # right-tailed
cv_upper <- qt(p = alpha, df = deg_freedom, lower.tail = F)
cv_legend <- round(cv_upper, 2)
cv_legend2 <- round(cv_upper * std_err + mu, 2)
pval <- round(pt(q = t, df = deg_freedom, lower.tail = F), 4)
} else { # two-tailed
cv_lower <- qt(p = alpha / 2, df = deg_freedom, lower.tail = T)
cv_upper <- qt(p = alpha / 2, df = deg_freedom, lower.tail = F)
cv_legend <- str_c('\u00B1', round(cv_upper, 2))
cv_legend2 <- str_c(round(cv_lower * std_err + mu, 2), ' & ', round(cv_upper * std_err + mu, 2))
pval_half <- round(pt(q = t, df = deg_freedom, lower.tail = t < 0), 4)
pval <- str_c(pval_half, ' + ', pval_half, ' = ', 2 * pval_half)
}
# Plot t-distribution
p <- ggplot(data.frame(x = -c(-4, 4)), aes(x)) +
ggtitle(plot_title) + xlab('') + ylab('') +
stat_function(fun = dt, args = list(df = deg_freedom), xlim = c(-4, 4))
# Shade rejection region using critical value
if (alternative != 'greater') {
p <- p + geom_vline(aes(xintercept = cv_lower, color = 'cval'),
linetype = 'dotted', size = 0.8, alpha = 0.8)
if (shade_rejection) {
p <- p + stat_function(fun = dt, args = list(df = deg_freedom),
xlim = c(-4, cv_lower),
geom = 'area', alpha = 0.3, fill = 'red')
}
if (shade_pval) {
p <- p + stat_function(fun = dt, args = list(df = deg_freedom),
xlim = c(-4, if_else(alternative == 'two.sided', -abs(t), t)),
geom = 'area', alpha = 0.3, fill = 'blue')
}
}
if (alternative != 'less') {
p <- p + geom_vline(aes(xintercept = cv_upper, color = 'cval'),
linetype = 'dotted', size = 0.8, alpha = 0.8)
if (shade_rejection) {
p <- p + stat_function(fun = dt, args = list(df = deg_freedom),
xlim = c(cv_upper, 4),
geom = 'area', alpha = 0.3, fill = 'red')
}
if (shade_pval) {
p <- p + stat_function(fun = dt, args = list(df = deg_freedom),
xlim = c(if_else(alternative == 'two.sided', abs(t), t), 4),
geom = 'area', alpha = 0.3, fill = 'blue')
}
}
# Legend text
legend_text <- c('t-statistics / p-value', 'critical value / alpha')
if (stacked) {
legend_text <- c(str_c('t-statistics: ', round(t, 2),
'\n(p-value: ', str_extract(pval, '[\\d.-]+$'), ')'),
str_c('Critical value: ', cv_legend,
'\n(alpha: ', round(alpha, 2), ')'))
}
stats_text <- c(str_c('t-statistics: ', round(t, 2)),
str_c('SE: ', round(std_err, 2)),
str_c('p-value: ', pval),
str_c('Critical value: ', cv_legend),
str_c('alpha: ', round(alpha, 2)))
if (!stacked) {
p <- p +
annotate('text', size = 9*5/14, x = 4.84, y = 0.14, hjust = 0,
label = 'bold(Statistics)', parse = T) +
annotate('text', size = 8*5/14, x = 4.89, y = 0:4 * -0.015 + 0.12, hjust = 0,
label = stats_text)
}
# Label plot
p <- p +
geom_vline(aes(xintercept = t, color = 'tstat'),
linetype = 'dotted', size = 0.8, alpha = 0.8) +
scale_x_continuous(sec.axis = sec_axis(trans = ~ . * std_err + mu)) +
scale_color_manual(name = if_else(stacked, 'Statistics', 'Legend'),
breaks = c('tstat', 'cval'),
labels = legend_text,
values = c(tstat = 'blue', cval = 'red')) +
theme(plot.title = element_text(size = 10, face = 'bold', hjust = 0.5),
plot.margin = unit(c(5.5, if_else(stacked, 5.5, 30), 5.5, 5.5), 'pt'),
legend.title = element_text(size = 9, face = 'bold'),
legend.text = element_text(size = 8)) +
coord_cartesian(xlim = c(-4, 4),
clip = 'off')
p
}2 Introduce IPEDS data
2.1 What/why IPEDS
What is IPEDS?
- ADD TEXT
Which IPEDS variables will we be using to teach statistical inference
- Variables about annual cost of attendance for full-time graduate programs
- average tuition (in-state vs. out-of-state)
- required fees (in-state vs. out-of-state)
- books and supplies; off-campus room and board; off-campus “other expenses”
- these are calculated for undergraduate students, but should be highly correlated with graduate students
- Note: none of these variables include the cost of healthcare
- Therefore, assuming you must pay for the cost of healthcare, the measures of cost of attendance (COA) we create will understate the total cost of graduate school + life
Why use IPEDS data rather than College Scorecard data? and why these variables?
- College Scorecard brings together variables from many data sources, including IPEDS
- The most interesting variables in College Scorecard are debt and earnings from degree programs, which come from Office of Federal Student Aid and IRS (I think)
- These university-level debt and earnings variables are measured as the mean (or median) of all students who graduated from the university
- so the underlying data (that we don’t see) are student-level; these student-level data are aggregated to the university-level to create measures mean (or median) debt/earnings at that university
- Using these measures would make teaching statistical inference more confusing (e.g., trying to test hypotheses about the population mean of mean earnings)
- Why teach inferential statistics using IPEDS measures about graduate tuition (and cost of attendance)?
- These measures are truly measured at the university-level
- Tuition and cost-of-attendance are the big drivers of student debt, so learning more about how tuition/cost-of-attendance varies across universities will help you make more informed decisions about pursuing graduate education in the future
2.2 Some definitions
Some definitions related to tuition, fees, expenses, etc; from the IPEDS “Student Charges for Full Academic Year” 2019-20 academic year data dictionary [LINK]:
- “Full-time student” (graduate)
- Graduate - A student enrolled for 9 or more semester credits, or 9 or more quarter credits, or students involved in thesis or dissertation preparation that are considered full time by the institution
- “Academic year”
- The period of time generally extending from September to June; usually equated to 2 semesters or trimesters, 3 quarters, or the period covered by a 4-1-4 plan
- “Tuition”
- Amount of money charged to students for instructional services. Tuition may be charged per term, per course, or per credit
- “In-state tuition”
- The tuition charged by institutions to those students who meet the state’s or institution’s residency requirements
- “Out-of-state tuition”
- The tuition charged by institutions to those students who do not meet the state’s or institution’s residency requirements
- “Required fees”
- Fixed sum charged to students for items not covered by tuition and required of such a large proportion of all students that the student who does NOT pay the charge is an exception
- “In-state fees”
- The fees charged by institutions to those students who meet the state’s or institution’s residency requirements.
- “Out-of-state fees”
- The fees charged by institutions to those students who do not meet the state’s or institution’s residency requirements
- “Books and supplies” (undergraduate)
- Do not include unusual costs for special groups of students (e.g., engineering or art majors), unless they constitute the majority of students at your institution
- “Room charges”
- The charges for an academic year for rooming accommodations for a typical student sharing a room with one other student.
- “Board charges”
- The charge for an academic year for meals, for a specified number of meals per week.
- “Other expenses”
- The amount of money (estimated by the financial aid office) needed by a student to cover expenses such as laundry, transportation, entertainment, and furnishings. (For the purpose of this survey room and board and tuition and fees are not included.)
Note: the IPEDS measures of full-time graduate tuition (both in-state and out-of-state) are “average” tuition price across different graduate degree programs (excluding “first-professional” degree programs like law and medicine)
- e.g., full-time annual tuition price for an MBA program is likely higher than that of an MA in education
- But in the context of using these measures to teach statistical inference, pretend these prices are not “averages” but rather the official tuition price
- make the same assumption for variables about fees, books/supplies, living expenses, etc.
Note: “required fees” do not include the cost of healthcare (I think)
2.3 Investigate df_ipeds_pop
Show variable labels
df_ipeds_pop %>% var_label()#> $instnm
#> [1] "Institution (entity) name"
#>
#> $unitid
#> [1] "Unique identification number of the institution"
#>
#> $opeid6
#> [1] "First 6 digits of OPEID"
#>
#> $opeid
#> [1] "Office of Postsecondary Education (OPE) ID Number"
#>
#> $control
#> [1] "Control of institution"
#>
#> $c15basic
#> [1] "Carnegie Classification 2015: Basic"
#>
#> $stabbr
#> [1] "State abbreviation"
#>
#> $city
#> [1] "City location of institution"
#>
#> $zip
#> [1] "ZIP code"
#>
#> $locale
#> [1] "Degree of urbanization (Urban-centric locale)"
#>
#> $region
#> [1] "Bureau of Economic Analysis (BEA) regions"
#>
#> $tuit_grad_res
#> [1] "In-state average tuition full-time graduates"
#>
#> $fee_grad_res
#> [1] "In-state required fees for full-time graduates"
#>
#> $tuit_grad_nres
#> [1] "Out-of-state average tuition full-time graduates"
#>
#> $fee_grad_nres
#> [1] "Out-of-state required fees for full-time graduates"
#>
#> $tuit_md_res
#> [1] "Medicine: In-state tuition"
#>
#> $fee_md_res
#> [1] "Medicine: In-state required fees"
#>
#> $tuit_md_nres
#> [1] "Medicine: Out-of-state tuition"
#>
#> $fee_md_nres
#> [1] "Medicine: Out-of-state required fees"
#>
#> $tuit_law_res
#> [1] "Law: In-state tuition"
#>
#> $fee_law_res
#> [1] "Law: In-state required fees"
#>
#> $tuit_law_nres
#> [1] "Law: Out-of-state tuition"
#>
#> $fee_law_nres
#> [1] "Law: Out-of-state required fees"
#>
#> $books_supplies
#> [1] "Books and supplies 2019-20"
#>
#> $roomboard_off
#> [1] "Off campus (not with family), room and board 2019-20"
#>
#> $oth_expense_off
#> [1] "Off campus (not with family), other expenses 2019-20"
#>
#> $tuitfee_grad_res
#> [1] "graduate, full-time, resident; avg tuition + required fees"
#>
#> $tuitfee_grad_nres
#> [1] "graduate, full-time, non-resident; avg tuition + required fees"
#>
#> $tuitfee_md_res
#> [1] "MD, full-time, state resident; avg tuition + required fees"
#>
#> $tuitfee_md_nres
#> [1] "MD, full-time, non-resident; avg tuition + required fees"
#>
#> $tuitfee_law_res
#> [1] "Law, full-time, state resident; avg tuition + required fees"
#>
#> $tuitfee_law_nres
#> [1] "Law, full-time, non-resident; avg tuition + required fees"
#>
#> $coa_grad_res
#> [1] "graduate, full-time, state resident COA; == tuition + fees + (ug) books/supplies + (ug) off-campus room and board + (ug) off-campus other expenses"
#>
#> $coa_grad_nres
#> [1] "graduate, full-time, non-resident COA; == tuition + fees + (ug) books/supplies + (ug) off-campus room and board + (ug) off-campus other expenses"
#>
#> $coa_md_res
#> [1] "MD, full-time, state resident COA; == tuition + fees + (ug) books/supplies + (ug) off-campus room and board + (ug) off-campus other expenses"
#>
#> $coa_md_nres
#> [1] "MD, full-time, non-resident COA; == tuition + fees + (ug) books/supplies + (ug) off-campus room and board + (ug) off-campus other expenses"
#>
#> $coa_law_res
#> [1] "Law, full-time, state resident COA; == tuition + fees + (ug) books/supplies + (ug) off-campus room and board + (ug) off-campus other expenses"
#>
#> $coa_law_nres
#> [1] "Law, full-time, non-resident COA; == tuition + fees + (ug) books/supplies + (ug) off-campus room and board + (ug) off-campus other expenses"Show value labels for variables that are labelled class (code note run)
df_ipeds_pop %>% select(control,locale,region,c15basic) %>% val_labels()Investigate data structure
- confirm one observation per unitid
df_ipeds_pop %>% group_by(unitid) %>% summarise(n_per_key=n()) %>% ungroup() %>% count(n_per_key)#> # A tibble: 1 x 2
#> n_per_key n
#> <int> <int>
#> 1 1 9912.3.1 Tuition/COA at UC campuses
Graduate, state residents
df_ipeds_pop %>%
# keep UC campuses
filter(unitid %in% c(110398,110635,110644,110653,110662,110671,110680,110699,110705,110714,445188,110699,110398)) %>%
select(instnm,unitid,tuit_grad_res,fee_grad_res,tuitfee_grad_res,books_supplies,roomboard_off,oth_expense_off,coa_grad_res)#> # A tibble: 9 x 9
#> instnm unitid tuit_grad_res fee_grad_res tuitfee_grad_res books_supplies roomboard_off oth_expense_off coa_grad_res
#> <chr> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
#> 1 University of California-Berkeley 110635 11442 2745 14187 870 14771 5359 35187
#> 2 University of California-Davis 110644 11442 2156 13598 1158 10588 4856 30200
#> 3 University of California-Irvine 110653 11442 1907 13349 1390 12861 5184 32784
#> 4 University of California-Los Angeles 110662 11442 1511 12953 1463 14303 5126 33845
#> 5 University of California-Riverside 110671 11442 2142 13584 1436 10986 4792 30798
#> 6 University of California-San Diego 110680 11442 2013 13455 1128 13681 4760 33024
#> 7 University of California-Santa Barbara 110705 11442 2112 13554 1184 12818 6045 33601
#> 8 University of California-Santa Cruz 110714 11442 2366 13808 1085 13216 5442 33551
#> 9 University of California-Merced 445188 11442 1765 13207 1016 8595 4909 27727Graduate, non-resident
df_ipeds_pop %>%
# keep UC campuses
filter(unitid %in% c(110398,110635,110644,110653,110662,110671,110680,110699,110705,110714,445188,110699,110398)) %>%
select(instnm,unitid,tuit_grad_nres,fee_grad_nres,tuitfee_grad_nres,books_supplies,roomboard_off,oth_expense_off,coa_grad_nres)#> # A tibble: 9 x 9
#> instnm unitid tuit_grad_nres fee_grad_nres tuitfee_grad_nres books_supplies roomboard_off oth_expense_off coa_grad_nres
#> <chr> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
#> 1 University of California-Berkeley 110635 26544 2745 29289 870 14771 5359 50289
#> 2 University of California-Davis 110644 26544 2156 28700 1158 10588 4856 45302
#> 3 University of California-Irvine 110653 26544 1907 28451 1390 12861 5184 47886
#> 4 University of California-Los Angeles 110662 26544 1511 28055 1463 14303 5126 48947
#> 5 University of California-Riverside 110671 26544 2142 28686 1436 10986 4792 45900
#> 6 University of California-San Diego 110680 26544 2013 28557 1128 13681 4760 48126
#> 7 University of California-Santa Barbara 110705 26544 2112 28656 1184 12818 6045 48703
#> 8 University of California-Santa Cruz 110714 26544 2366 28910 1085 13216 5442 48653
#> 9 University of California-Merced 445188 26544 1765 28309 1016 8595 4909 42829MD, state resident
df_ipeds_pop %>%
# keep UC campuses
filter(unitid %in% c(110398,110635,110644,110653,110662,110671,110680,110699,110705,110714,445188,110699,110398)) %>%
select(instnm,unitid,tuit_md_res,fee_md_res,tuitfee_md_res,books_supplies,roomboard_off,oth_expense_off,coa_md_res)#> # A tibble: 9 x 9
#> instnm unitid tuit_md_res fee_md_res tuitfee_md_res books_supplies roomboard_off oth_expense_off coa_md_res
#> <chr> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
#> 1 University of California-Berkeley 110635 34521 2745 37266 870 14771 5359 58266
#> 2 University of California-Davis 110644 34526 2156 36682 1158 10588 4856 53284
#> 3 University of California-Irvine 110653 34527 1907 36434 1390 12861 5184 55869
#> 4 University of California-Los Angeles 110662 36374 1511 37885 1463 14303 5126 58777
#> 5 University of California-Riverside 110671 35430 2142 37572 1436 10986 4792 54786
#> 6 University of California-San Diego 110680 34977 2013 36990 1128 13681 4760 56559
#> 7 University of California-Santa Barbara 110705 NA NA NA 1184 12818 6045 NA
#> 8 University of California-Santa Cruz 110714 NA NA NA 1085 13216 5442 NA
#> 9 University of California-Merced 445188 NA NA NA 1016 8595 4909 NAMD, non-resident
df_ipeds_pop %>%
# keep UC campuses
filter(unitid %in% c(110398,110635,110644,110653,110662,110671,110680,110699,110705,110714,445188,110699,110398)) %>%
select(instnm,unitid,tuit_md_nres,fee_md_nres,tuitfee_md_nres,books_supplies,roomboard_off,oth_expense_off,coa_md_nres)#> # A tibble: 9 x 9
#> instnm unitid tuit_md_nres fee_md_nres tuitfee_md_nres books_supplies roomboard_off oth_expense_off coa_md_nres
#> <chr> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
#> 1 University of California-Berkeley 110635 46766 2745 49511 870 14771 5359 70511
#> 2 University of California-Davis 110644 46771 2156 48927 1158 10588 4856 65529
#> 3 University of California-Irvine 110653 46772 1907 48679 1390 12861 5184 68114
#> 4 University of California-Los Angeles 110662 48619 1511 50130 1463 14303 5126 71022
#> 5 University of California-Riverside 110671 47675 2142 49817 1436 10986 4792 67031
#> 6 University of California-San Diego 110680 47222 2013 49235 1128 13681 4760 68804
#> 7 University of California-Santa Barbara 110705 NA NA NA 1184 12818 6045 NA
#> 8 University of California-Santa Cruz 110714 NA NA NA 1085 13216 5442 NA
#> 9 University of California-Merced 445188 NA NA NA 1016 8595 4909 NALaw, state resident
df_ipeds_pop %>%
# keep UC campuses
filter(unitid %in% c(110398,110635,110644,110653,110662,110671,110680,110699,110705,110714,445188,110699,110398)) %>%
select(instnm,unitid,tuit_law_res,fee_law_res,tuitfee_law_res,books_supplies,roomboard_off,oth_expense_off,coa_law_res)#> # A tibble: 9 x 9
#> instnm unitid tuit_law_res fee_law_res tuitfee_law_res books_supplies roomboard_off oth_expense_off coa_law_res
#> <chr> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
#> 1 University of California-Berkeley 110635 49242 2745 51987 870 14771 5359 72987
#> 2 University of California-Davis 110644 47332 2156 49488 1158 10588 4856 66090
#> 3 University of California-Irvine 110653 45578 1907 47485 1390 12861 5184 66920
#> 4 University of California-Los Angeles 110662 45578 1511 47089 1463 14303 5126 67981
#> 5 University of California-Riverside 110671 NA NA NA 1436 10986 4792 NA
#> 6 University of California-San Diego 110680 NA NA NA 1128 13681 4760 NA
#> 7 University of California-Santa Barbara 110705 NA NA NA 1184 12818 6045 NA
#> 8 University of California-Santa Cruz 110714 NA NA NA 1085 13216 5442 NA
#> 9 University of California-Merced 445188 NA NA NA 1016 8595 4909 NALaw, non-resident
df_ipeds_pop %>%
# keep UC campuses
filter(unitid %in% c(110398,110635,110644,110653,110662,110671,110680,110699,110705,110714,445188,110699,110398)) %>%
select(instnm,unitid,tuit_law_nres,fee_law_nres,tuitfee_law_nres,books_supplies,roomboard_off,oth_expense_off,coa_law_nres)#> # A tibble: 9 x 9
#> instnm unitid tuit_law_nres fee_law_nres tuitfee_law_nres books_supplies roomboard_off oth_expense_off coa_law_nres
#> <chr> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
#> 1 University of California-Berkeley 110635 52571 2745 55316 870 14771 5359 76316
#> 2 University of California-Davis 110644 56435 2156 58591 1158 10588 4856 75193
#> 3 University of California-Irvine 110653 51641 1907 53548 1390 12861 5184 72983
#> 4 University of California-Los Angeles 110662 51641 1511 53152 1463 14303 5126 74044
#> 5 University of California-Riverside 110671 NA NA NA 1436 10986 4792 NA
#> 6 University of California-San Diego 110680 NA NA NA 1128 13681 4760 NA
#> 7 University of California-Santa Barbara 110705 NA NA NA 1184 12818 6045 NA
#> 8 University of California-Santa Cruz 110714 NA NA NA 1085 13216 5442 NA
#> 9 University of California-Merced 445188 NA NA NA 1016 8595 4909 NA2.3.2 Tuition/COA at fancy privates
Some fancy-pants private universities
- USC
unitid == 123961
- Stanford
unitid == 243744
- Columbia
unitid == 190150
- Columbia, Teacher’s College
unitid == 196468
- NYU
unitid == 193900
- Harvard
unitid == 166027
- Vanderbilt
unitid == 221999
- University of Pennsylvania
unitid == 215062
- Northwestern University
unitid == 147767
- Johns Hopkins University
unitid == 162928
Graduate students
- Note: no difference in prices between resident and non-resident
- Note: Some universities did not provide
books_supplies,roomboard_off,oth_expense_off - Some ludicrously unrealistic values for
roomboard_offandoth_expense_off
# In-state
df_ipeds_pop %>%
# keep private fancy pants
filter(unitid %in% c(123961,243744,190150,196468,193900,166027,221999,215062,147767,162928)) %>%
select(instnm,unitid,tuit_grad_res,fee_grad_res,tuitfee_grad_res,books_supplies,roomboard_off,oth_expense_off,coa_grad_res)#> # A tibble: 6 x 9
#> instnm unitid tuit_grad_res fee_grad_res tuitfee_grad_res books_supplies roomboard_off oth_expense_off coa_grad_res
#> <chr> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
#> 1 University of Southern California 123961 46272 835 47107 1200 15912 2152 66371
#> 2 Northwestern University 147767 56067 500 56567 1638 17019 3306 78530
#> 3 Johns Hopkins University 162928 55350 2206 57556 1250 12235 1118 72159
#> 4 Columbia University in the City of New York 190150 47600 2368 49968 1294 17955 5220 74437
#> 5 New York University 193900 33408 2154 35562 752 18684 2758 57756
#> 6 University of Pennsylvania 215062 36254 3928 40182 1358 14406 1946 57892# Out-of-state
df_ipeds_pop %>%
# keep private fancy pants
filter(unitid %in% c(123961,243744,190150,196468,193900,166027,221999,215062,147767,162928)) %>%
select(instnm,unitid,tuit_grad_nres,fee_grad_nres,tuitfee_grad_nres,books_supplies,roomboard_off,oth_expense_off,coa_grad_nres)#> # A tibble: 6 x 9
#> instnm unitid tuit_grad_nres fee_grad_nres tuitfee_grad_nres books_supplies roomboard_off oth_expense_off coa_grad_nres
#> <chr> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
#> 1 University of Southern California 123961 46272 835 47107 1200 15912 2152 66371
#> 2 Northwestern University 147767 56067 500 56567 1638 17019 3306 78530
#> 3 Johns Hopkins University 162928 55350 2206 57556 1250 12235 1118 72159
#> 4 Columbia University in the City of New York 190150 47600 2368 49968 1294 17955 5220 74437
#> 5 New York University 193900 33408 2154 35562 752 18684 2758 57756
#> 6 University of Pennsylvania 215062 36254 3928 40182 1358 14406 1946 57892MD students
df_ipeds_pop %>%
# keep private fancy pants
filter(unitid %in% c(123961,243744,190150,196468,193900,166027,221999,215062,147767,162928)) %>%
select(instnm,unitid,tuit_md_res,fee_md_res,tuitfee_md_res,books_supplies,roomboard_off,oth_expense_off,coa_md_res)#> # A tibble: 6 x 9
#> instnm unitid tuit_md_res fee_md_res tuitfee_md_res books_supplies roomboard_off oth_expense_off coa_md_res
#> <chr> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
#> 1 University of Southern California 123961 64538 3065 67603 1200 15912 2152 86867
#> 2 Northwestern University 147767 62088 828 62916 1638 17019 3306 84879
#> 3 Johns Hopkins University 162928 54900 4909 59809 1250 12235 1118 74412
#> 4 Columbia University in the City of New York 190150 62980 2340 65320 1294 17955 5220 89789
#> 5 New York University 193900 55645 3495 59140 752 18684 2758 81334
#> 6 University of Pennsylvania 215062 59910 5378 65288 1358 14406 1946 82998Law students
df_ipeds_pop %>%
# keep private fancy pants
filter(unitid %in% c(123961,243744,190150,196468,193900,166027,221999,215062,147767,162928)) %>%
select(instnm,unitid,tuit_law_res,fee_law_res,tuitfee_law_res,books_supplies,roomboard_off,oth_expense_off,coa_law_res)#> # A tibble: 6 x 9
#> instnm unitid tuit_law_res fee_law_res tuitfee_law_res books_supplies roomboard_off oth_expense_off coa_law_res
#> <chr> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
#> 1 University of Southern California 123961 66306 885 67191 1200 15912 2152 86455
#> 2 Northwestern University 147767 66506 300 66806 1638 17019 3306 88769
#> 3 Johns Hopkins University 162928 NA NA NA 1250 12235 1118 NA
#> 4 Columbia University in the City of New York 190150 69896 2464 72360 1294 17955 5220 96829
#> 5 New York University 193900 66036 2676 68712 752 18684 2758 90906
#> 6 University of Pennsylvania 215062 63610 4996 68606 1358 14406 1946 863162.3.3 Some basic descriptive statistics
Number of institutions by type
# Number of institutions
df_ipeds_pop %>% count()#> # A tibble: 1 x 1
#> n
#> <int>
#> 1 991# Number of institutions by Carnegie type
df_ipeds_pop %>% count(unclass(c15basic)) %>% as_factor()#> # A tibble: 6 x 2
#> `unclass(c15basic)` n
#> <dbl> <int>
#> 1 15 102
#> 2 16 98
#> 3 17 95
#> 4 18 374
#> 5 19 199
#> 6 20 123# Number of institutions by public/private
df_ipeds_pop %>% count(unclass(control)) %>% as_factor()#> # A tibble: 3 x 2
#> `unclass(control)` n
#> <dbl> <int>
#> 1 1 459
#> 2 2 487
#> 3 3 45# number of institutions by public/private and carnegie type
df_ipeds_pop %>% count(unclass(control),unclass(c15basic)) %>% as_factor()#> # A tibble: 16 x 3
#> `unclass(control)` `unclass(c15basic)` n
#> <dbl> <dbl> <int>
#> 1 1 15 79
#> 2 1 16 74
#> 3 1 17 38
#> 4 1 18 159
#> 5 1 19 68
#> 6 1 20 41
#> 7 2 15 23
#> 8 2 16 24
#> 9 2 17 50
#> 10 2 18 193
#> 11 2 19 121
#> 12 2 20 76
#> 13 3 17 7
#> 14 3 18 22
#> 15 3 19 10
#> 16 3 20 6# number of institutions by level of urbanization
df_ipeds_pop %>% count(unclass(locale)) %>% as_factor()#> # A tibble: 12 x 2
#> `unclass(locale)` n
#> <dbl> <int>
#> 1 11 252
#> 2 12 142
#> 3 13 147
#> 4 21 195
#> 5 22 25
#> 6 23 27
#> 7 31 25
#> 8 32 84
#> 9 33 65
#> 10 41 17
#> 11 42 8
#> 12 43 4# number of institutions by public/private and level of urbanization
df_ipeds_pop %>% count(unclass(control),unclass(locale)) %>% as_factor() %>% head(10)#> # A tibble: 10 x 3
#> `unclass(control)` `unclass(locale)` n
#> <dbl> <dbl> <int>
#> 1 1 11 95
#> 2 1 12 67
#> 3 1 13 78
#> 4 1 21 57
#> 5 1 22 17
#> 6 1 23 16
#> 7 1 31 13
#> 8 1 32 51
#> 9 1 33 54
#> 10 1 41 6Tuition+fees by public/private and Carnegie type
df_ipeds_pop %>% group_by(unclass(control),unclass(c15basic)) %>%
summarize(
sample_size = n(),
n_nonmiss_tuitfee_res = sum(!is.na(tuitfee_grad_res)),
mean_tuitfee_res = mean(tuitfee_grad_res, na.rm = TRUE),
n_nonmiss_tuitfee_nres = sum(!is.na(tuitfee_grad_nres)),
mean_tuitfee_nres = mean(tuitfee_grad_nres, na.rm = TRUE),
) %>% as_factor()#> # A tibble: 16 x 7
#> # Groups: unclass(control) [3]
#> `unclass(control)` `unclass(c15basic)` sample_size n_nonmiss_tuitfee_res mean_tuitfee_res n_nonmiss_tuitfee_nres mean_tuitfee_nres
#> <dbl> <dbl> <int> <int> <dbl> <int> <dbl>
#> 1 1 15 79 79 13635. 79 28419.
#> 2 1 16 74 74 11010. 74 22425.
#> 3 1 17 38 38 9735. 38 19672.
#> 4 1 18 159 159 10241. 159 18280.
#> 5 1 19 68 68 8909. 68 16047.
#> 6 1 20 41 41 9567. 41 17019.
#> 7 2 15 23 23 45740. 23 45740.
#> 8 2 16 24 24 30842. 24 30842.
#> 9 2 17 50 50 20525. 50 20525.
#> 10 2 18 193 193 16842. 193 16842.
#> 11 2 19 121 121 14498. 121 14498.
#> 12 2 20 76 76 14159. 76 14159.
#> 13 3 17 7 7 13066. 7 13066.
#> 14 3 18 22 22 14110. 22 14110.
#> 15 3 19 10 10 16398. 10 16398.
#> 16 3 20 6 6 19569. 6 19933.Cost of attendance by public/private and Carnegie type
df_ipeds_pop %>% group_by(unclass(control),unclass(c15basic)) %>%
summarize(
sample_size = n(),
n_nonmiss_coa_res = sum(!is.na(coa_grad_res)),
mean_coa_res = mean(coa_grad_res, na.rm = TRUE),
n_nonmiss_coa_nres = sum(!is.na(coa_grad_nres)),
mean_coa_nres = mean(coa_grad_nres, na.rm = TRUE),
) %>% as_factor()#> # A tibble: 16 x 7
#> # Groups: unclass(control) [3]
#> `unclass(control)` `unclass(c15basic)` sample_size n_nonmiss_coa_res mean_coa_res n_nonmiss_coa_nres mean_coa_nres
#> <dbl> <dbl> <int> <int> <dbl> <int> <dbl>
#> 1 1 15 79 79 30388. 79 45173.
#> 2 1 16 74 74 27004. 74 38419.
#> 3 1 17 38 38 25972. 38 35909.
#> 4 1 18 159 159 26449. 159 34487.
#> 5 1 19 68 68 24017. 68 31155.
#> 6 1 20 41 41 24417. 41 31870.
#> 7 2 15 23 23 64642. 23 64642.
#> 8 2 16 24 24 48278. 24 48278.
#> 9 2 17 50 50 37475. 50 37475.
#> 10 2 18 193 193 32363. 193 32363.
#> 11 2 19 121 121 29331. 121 29331.
#> 12 2 20 76 76 28431. 76 28431.
#> 13 3 17 7 7 25133 7 25133
#> 14 3 18 22 22 28743. 22 28743.
#> 15 3 19 10 10 32577. 10 32577.
#> 16 3 20 6 6 35184. 6 35548.3 Distributions
Definitions
- Probability distribution
In probability theory and statistics, a probability distribution is the mathematical function that gives the probabilities of occurrence of different possible outcomes
A little easier to get your head around “frequency distribution” and “relative frequency distribution”
- Frequency distribution
- for each value of a variable, the number of observations that have that value
- Relative frequency distribution
- for each value of a variable, the percent of observations that have that value (the number of observations with that value divided by the total number of observations)
We can show frequency distributions (or relative frequency distributions) as a table or as a graph
- Frequency distribution of variable
tuit_grad_res(for first few values)
df_ipeds_pop %>% count(tuitfee_grad_res) %>% head(10)#> # A tibble: 10 x 2
#> tuitfee_grad_res n
#> <dbl> <int>
#> 1 945 1
#> 2 2320 1
#> 3 3552 3
#> 4 3620 1
#> 5 3710 1
#> 6 3915 1
#> 7 4140 1
#> 8 4146 1
#> 9 4172 1
#> 10 4356 1We can also visualize the distribution of the variable tuitfee_grad_res.
- Instead of “visualize,” sometimes we will say “plot” or “draw”
Below, we call a function named plot_distribution() that we will call on below to plot variables. Our plot_distribution() function calls functions from the ggplot2 library, which is part of the tidyverse
- (Props to Crystal Han (data scientist extraordinaire) for writing this function – and some other functions I use below – in short order as I was way behind in writing this lecture!)
Call the plot_distribution() function to plot the variable tuitfee_grad_res from the data frame df_ipeds_pop. Essentially, the plot_distribution() function creates the following things:
- draws a “histogram” of the variable
- draws a “kernel density estimate,” which – according to the help file for the
geom_density()function – is a “smoothed version of the histogram. This is a useful alternative to the histogram for continuous data that comes from an underlying smooth distribution” - prints the value of the following descriptive statistics: mean, median, standard deviation
- draws dotted vertical lines to indicate the value of the mean and the median
plot_distribution(df_ipeds_pop, 'tuitfee_grad_res')
3.1 Normal, right-skew, left-skew
The distributions of variables fall into a few general categories
- we will talk about normally distributed variables, right-skewed variables, and left-skewed variables
In order to explain these concepts, we will generate a normally distributed variable, a right-skewed variable, and a left-skewed variable
- We will pretend that each variable contains all observations in the population for some population of interest
- Usually, we don’t know the population; we collect a sample from that population and use the sample data to make statements (inferences) about the unknown population
- we will put each of these variables in a dataset called
def_generated_pop - don’t worry about below code
We can refer to variables in the data frame df_generated_pop using the following syntax:
data_frame_name$variable_name- e.g.,
df_generated_pop$norm_dist
Examine the variable norm_dist in the data frame df_generated_pop
length(df_generated_pop$norm_dist) # length() = number of "elements" = number of observations#> [1] 10000#mean
mean(df_generated_pop$norm_dist, na.rm = TRUE)#> [1] 49.98631Normal distribution
- a normal distribution is “bell shaped”; has symmetric “tails”
We generated a variable df_generated_pop$norm_dist that has a normal distribution and then plot the variable to visualize what a normal distribution looks
- Descriptive statistics about the variable
df_generated_pop$norm_dist- It has a mean of 49.99
- It has a standard deviation of 4.99
- Standard deviation is a measure of how far away from the mean observations tend to be
- we can interpet this standard deviation as follows: on average, observations are 4.99 away from the mean of 49.99
We can also visualize the variable df_generated_pop$norm_dist, as shown below. Note the following:
- symmetric, “bell” shape
- the mean is (nearly) identical to the median
plot_distribution(df_generated_pop, 'norm_dist')
Skewed distribution
- for a variable with a “skewed distribution” one “tail” is longer than the other
- “outliers” are often the cause of skewed distributions
- outliers are observations that are much higher or much lower than most observations
- skewed distributions are either “left skewed” or “right skewed”
Right skewed distributions
- The right “tail” is longer than the left due to the presence of positive outliers, defined as observations with very high values compared to most observations
- so there are more positive outliers than you would expect in a bell (normal) shaped variable
- these positive outliers increase the value of the mean, such that the value of the mean is higher than the value of the median
- mean > median
- real-world variables that tend to be right-skewed
- income; enrollment size, city population
We generated the right-skewed variable df_generated_pop$rskew_dist
- This variable has the following descriptive statistics
- mean equals 0.284
- median equals 0.264
- standard deviation equals 0.159
- interpretation: on average, observations are 0.159 away from the mean of 0.284
- We can also visualize the right-skwed variable
df_generated_pop$rskew_dist, as shown below. Note the following:- The right “tail” is longer than the left tail
- the mean is larger than the median
plot_distribution(df_generated_pop, 'rskew_dist')
Left skewed distributions
- The left tail is longer than right tail, usually due to the presence of more negative outliers than would be expected in a bell shaped variable
- where negative outliers are defined as observations with very low values (e.g., extreme negative values) compared to most observations
- these negative outliers decrease the value of the mean, such that the value of the mean is lower than the value of the median
- In social science research left-skewed variables are less common than right-skewed variables
We generated the left-skewed variable df_generated_pop$rskew_dist
- This variable has the following descriptive statistics
- mean equals 0.716
- median equals 0.736
- standard deviation equals 0.159
- interpretation: on average, observations are 0.159 away from the mean of 0.716
- We can also visualize the left-skwed variable
df_generated_pop$rskew_dist, as shown below. Note the following:- The left “tail” is longer than the right tail
- the mean is smaller than the median
Create and plot left-skewed variable
plot_distribution(df_generated_pop, 'lskew_dist')
Having plotted generated variables that have normal, right skewed, and left skewed distributions, respectively, let’s plot the real-life variable tuitfee_grad_res from the data frame df_ipeds_pop (below). How would you diagnose the distribution of tuitfee_grad_res?
plot_distribution(df_ipeds_pop, 'tuitfee_grad_res')
How would you diagnose the distribution of tuitfee_grad_res?
- right tail longer than the left
- mean of 14602 is larger than the median of 12343
- so this is a right skewed variable
3.2 Properties of the normal distribution
The normal distribution – symmetric, bell-shaped, mean equal to the median - has very useful properties that make it the basis of inferential statistics
- if we can reasonably assume a variable has a normal distribution, then we know a lot about that variable
plot_distribution(df_generated_pop, 'norm_dist')
3.2.1 Normal distribution and “the empirical rule”
Recall our primary measure of dispersion, standard deviation, which measures how far away individual observations tend to be from the mean.
- Here, for some variable \(x\), I’ll give the formula for sample standard deviation, \(\hat\sigma_x\), as opposed to the formula for population standard deviation, \(\sigma_x\)
- sometimes people will also refer to sample standard deviation of \(x\) as \(s_x\)
- \(\hat\sigma_x = \sqrt{\frac{\sum_{i=1}^n (x_i - \overline{x})^2}{n-1}}\)
The empirical rule states that if variable has an approximately normal distribution (i.e., approximately “bell shaped”) then:
- About 68% of obs fall within one std. dev of mean
- i.e., between \(x - \hat{\sigma{x}}\) and \(x + \hat{\sigma{x}}\)
- About 95% of obs fall within two std. dev of mean
- i.e., between \(x - 2\hat\sigma_x\) and \(x + 2\hat\sigma_x\)
- About 99% of obs fall within three std. dev of mean
- i.e., between \(x - 3\hat\sigma_x\) and \(x + 3\hat\sigma_x\)
Why is the empirical rule so important for inferential statistics?
- basically, if a variable has an approximately normal distribution, then we know how likely it would be to observe a variable that is a certain number of standard deviations away from the mean
- for example:
- only about 2.5% of observations have a value higher than two standard deviations or more from the mean;
- the variable
norm_disthas a mean of about50and a standard deviation of about5, so the value of40would be about two standard deviations below the mean. the empirical rule tells us that only about 2.5% of observations would have a value less than40
- you might say, but most real-life variables are unlikely to have a normal distribution
- True! But the “sampling distribution” – discussed below – which is the basis for all inferential statistics/hypothesis testing, always has a normal distribution so long as our sample size is large enough
3.2.2 Z-scores
The “z-score” of an observation is the number of standard deviations away from the mean
- can define a z-score in terms of a population or a sample from that population; we’ll introduce z-score in the context of sample data
z-score formula
- where \(x\) is some variable of interest; subscript \(i\) refers to observations
- \(z_i = (x_i - \bar{x})/(\hat{\sigma}_x)\)
- in words:
- z score for observation \(i\) equals the difference between the observation \(x_i\) and the mean \(\bar{x}\) divided by the standard deviation \(\hat{\sigma}_x\)
- Intuition behing z-score
- it’s just the difference between an observation value and the mean value, scaled in terms of standard deviations
- that’s why we say that the z-score represents the number of standard deviations away from the mean
Calculating z-score for the variable norm_dist from data frame df_generated_pop
# components of z-score
mean(df_generated_pop$norm_dist, na.rm = TRUE)#> [1] 49.98631sd(df_generated_pop$norm_dist, na.rm = TRUE)#> [1] 4.991961#create new variable z_norm_dist
df_generated_pop <- df_generated_pop %>% mutate(
z_norm_dist = (norm_dist - mean(norm_dist, na.rm = TRUE))/sd(norm_dist, na.rm = TRUE)
)
#list a few observations
df_generated_pop %>% select(norm_dist,z_norm_dist) %>% head(10)#> # A tibble: 10 x 2
#> norm_dist z_norm_dist
#> <dbl> <dbl>
#> 1 43.1 -1.38
#> 2 50.2 0.0411
#> 3 46.2 -0.762
#> 4 51.1 0.215
#> 5 57.1 1.43
#> 6 53.7 0.748
#> 7 53.5 0.704
#> 8 48.9 -0.227
#> 9 51.0 0.200
#> 10 56.0 1.21# mean of z-score variable
round(mean(df_generated_pop$z_norm_dist, na.rm = TRUE), digits = 4)#> [1] 0Plot the new z-score variable, which has:
- mean of about
0 - standard deviation of about
1
plot_distribution(df_generated_pop, 'z_norm_dist')
Recall empirical rule for normally distributed variables
- About 68% of obs fall within one std. dev of mean
- i.e., between \(x - z\) and \(x + z\)
- About 95% of obs fall within two std. dev of mean
- i.e., between \(x - 2z\) and \(x + 2z\)
- About 99% of obs fall within three std. dev of mean
- i.e., between \(x - 3z\) and \(x + 3z\)
Show picture from old “chicken scratch” pdf conveying ideas about z-score and empirical rule
Delete variable z_norm_dist
df_generated_pop$z_norm_dist <- NULL3.3 Standard normal distribution
Definition:
- a bell-shaped (i.e., normal) distribution that has a mean of
0and a standard deviation of1
Above, we created a variable stdnorm_dist in the data frame df_generated_pop that has a standard normal distribution. Let’s investigate and plot this variable:
mean(df_generated_pop$stdnorm_dist, na.rm = TRUE)#> [1] -0.002737966sd(df_generated_pop$stdnorm_dist, na.rm = TRUE)#> [1] 0.9983922plot_distribution(df_generated_pop, 'stdnorm_dist')
Commentary about standard normal distribution:
- The value of each observation is already in terms of z-scores
- That is, the value of each observation shows how many standard deviations it is from the mean
- Question: if the variable has a standard normal distribution, would it be likely to see an observation with a value of 3?
- No. because a value of 3 would mean that the observation is three standard deviations greater than the mean. we know that for any variable with a normal distribution, less than 1% of observations have a value that is three standard deviations greater than the mean.
Question: If we have a variable with a roughly normal distribution (i.e., symmetrical), how could we transform it into a variable with a standard normal distribution (i.e., symmetrical with mean=0 and std deviation=1)?
- Answer: create new variable that is the z-score associated with each value
- This is exactly, what we did above (shown again here)
# create variable that has a standard normal distribution from a variable that has a normal distribution
df_generated_pop %>% mutate(
z_norm_dist = (norm_dist - mean(norm_dist, na.rm = TRUE))/sd(norm_dist, na.rm = TRUE)
) %>% head(10)#> # A tibble: 10 x 5
#> norm_dist rskew_dist lskew_dist stdnorm_dist z_norm_dist
#> <dbl> <dbl> <dbl> <dbl> <dbl>
#> 1 43.1 0.0827 0.917 -1.39 -1.38
#> 2 50.2 0.294 0.706 0.0383 0.0411
#> 3 46.2 0.156 0.844 -0.763 -0.762
#> 4 51.1 0.331 0.669 0.212 0.215
#> 5 57.1 0.651 0.349 1.43 1.43
#> 6 53.7 0.460 0.540 0.744 0.748
#> 7 53.5 0.448 0.552 0.700 0.704
#> 8 48.9 0.242 0.758 -0.229 -0.227
#> 9 51.0 0.327 0.673 0.197 0.200
#> 10 56.0 0.588 0.412 1.21 1.214 Sampling distribution
The “sampling distribution” is the fundamental concept of inferential statistics
Briefly, recall the goal of inferential statistics:
- We want to make statements about population parameters, for example the population mean of some variable \(x\), \(\mu_x\)
- But we usually cannot obtain data on the entire population
- Therefore, we collect a representative (random) sample from the population
- We calculate “estimates” based on this sample data. These estimates are our best guess abut the value of population parameters
- For example, the sample mean \(\bar{x}\) is our best guess of the population mean \(\mu_x\)
Usually, we collect a single sample from the population. How do we know if the sample we collected is representative of the underlying population we want to make statements about?
- This is a problem that statisticians have thought a lot about
For example, for our variable norm_dist from the data frame df_generated_pop, we randomly draw 30 observations from a population of 10,000 observations
set.seed(321)
norm_dist_s1 <- sample(x = df_generated_pop$norm_dist, size = 30)
mean(df_generated_pop$norm_dist)#> [1] 49.98631mean(norm_dist_s1) # mean of our sample#> [1] 50.19565But, what if we had obtained a different random sample?
set.seed(123)
norm_dist_s1 <- sample(x = df_generated_pop$norm_dist, size = 30)
mean(df_generated_pop$norm_dist)#> [1] 49.98631mean(norm_dist_s1) # mean of our sample#> [1] 49.22062# remove object norm_dist_s1
rm(norm_dist_s1)So we can see that the sample mean, \(\bar{x}\), changes from sample to sample
Imagine if we take 1,000 random samples of size n (e.g., 30) from a population
- for each random sample, we calculate the sample mean, and record the value of the sample mean
- we would have 1,000 observations, where each observation is the value of a sample mean
- If we plotted these 1,000 observations, it would give us a distribution of sample means
- More specifically, this would give us the “sampling distribution” of the sample mean
Sampling distribution (of the sample mean)
- The sampling distribution of the sample mean is a relative frequency distribution where each observation is the sample mean of a single random sample from a population
- The sampling distribution shows how the value of the sample mean varies from sample to sample
Excellent website for understanding how the sampling distribution works
- This very useful website does interactive simulations that show how the sampling distribution works LINK
- Please spend 5 minutes playing around with this website; this is really the most important concept in statistics
The sampling distribution of any statistic
- So far, we have discussed the sampling distribution of the sample mean, but a sampling distribution can be created for any sample statistic (e.g., median, min, max, regression coefficient)
- Once we get to the unit on regression, we’ll be thinking about the sampling distribution of a regression coefficient. But the underlying concepts will be exactly the same as the sampling distribution of the sample mean
4.1 Plot sampling distribution
Plotting a single random sample and a sampling distribution
- We are pretending our two data frames
df_generated_popanddf_ipeds_popcontain all observations in the population - For each of these datasets we will create a version that is a single random sample
- We will plot this single random sample
- Then we will take a large number of random samples, calculate the sample mean for each, and plot these sample means to create the sampling distribution
Plot the variable norm_dist from the sample dataset
plot_distribution(df_generated_sample, 'norm_dist')
Plot the variable tuitfee_grad_res from the sample IPEDS dataset
plot_distribution(df_ipeds_sample, 'tuitfee_grad_res')
Plot the sampling distribution of of the sample mean for the normally distributed variable norm_dist
plot_distribution(df_generated_pop, 'norm_dist', sampling_dist = T)
Plot the sampling distribution of of the sample mean for the right-skewed variable rskew_dist
plot_distribution(df_generated_pop, 'rskew_dist', sampling_dist = T)
Create visualization that stacks three plots on top of one another: 1= population distribution; 2 = distribution of a single random sample; 3 = sampling distribution
- do this for the normally distributed variable
norm_dist
plot_distribution(df_generated_pop, 'norm_dist', plot_title = 'Population distribution') +
plot_distribution(df_generated_sample, 'norm_dist', plot_title = 'Single sample distribution') +
plot_distribution(df_generated_pop, 'norm_dist', sampling_dist = T,
plot_title = 'Sampling distribution') +
plot_layout(ncol = 1)
Create visualization that stacks three plots on top of one another: 1= population distribution; 2 = distribution of a single random sample; 3 = sampling distribution
- do this for the right-skewed variable
rskew_dist - Notice that even though the underlying population is right skewed and the sample is right-skewed, the sampling distribution has a normal distribution
plot_distribution(df_generated_pop, 'rskew_dist', plot_title = 'Population distribution') +
plot_distribution(df_generated_sample, 'rskew_dist', plot_title = 'Single sample distribution') +
plot_distribution(df_generated_pop, 'rskew_dist', sampling_dist = T,
plot_title = 'Sampling distribution') +
plot_layout(ncol = 1)
Create visualization that stacks three plots on top of one another: 1= population distribution; 2 = distribution of a single random sample; 3 = sampling distribution
- do this for the right-skewed variable
tuitfee_grad_resfrom IPEDS - Notice that even though the underlying population is right skewed and the sample is right-skewed, the sampling distribution has a normal distribution
plot_distribution(df_ipeds_pop, 'tuitfee_grad_res', plot_title = 'Population distribution') +
plot_distribution(df_ipeds_sample, 'tuitfee_grad_res', plot_title = 'Single sample distribution') +
plot_distribution(df_ipeds_pop, 'tuitfee_grad_res', sampling_dist = T,
plot_title = 'Sampling distribution') +
plot_layout(ncol = 1)
4.2 Central limit theorem
Central limit theorem
- For random sampling with a large sample size, \(n\), the sampling distribution of the sample mean is approximately normally distributed
- In other words, no matter what the distribution of the variable, the sampling distribution will have a normal distribution
What counts as a “large” sample size?
- some say
n=30or more
Why is central limit theorem important?
- We conducting hypothesis tests about a population parameter (e.g., about a population mean, about a population regression coefficient), based on the sampling distribution of the relevant sample statistic
- even if the underlying variable of interest has a skewed population, the sampling distribution of the sample statistic (e.g., sample mean) will have a normal distribution
- if the sampling distribution has a normal distribution, then we know the percent of observations that are a certain number of standard deviations from the mean
Show central limit theorem using interactive simulation
Show central limit theorem using a very skewed variable: non-resident, grad school cost of attendance
plot_distribution(df_ipeds_pop, 'tuitfee_grad_nres', plot_title = 'Population distribution') +
plot_distribution(df_ipeds_sample, 'tuitfee_grad_nres', plot_title = 'Single sample distribution') +
plot_distribution(df_ipeds_pop, 'tuitfee_grad_nres', sampling_dist = T,
plot_title = 'Sampling distribution') +
plot_layout(ncol = 1)
4.3 Standard error
Sample standard deviation, of some variable \(Y\)
- Sample standard deviation, \(\hat{\sigma}_Y\), is the average distance between a random observation and the sample mean \(\bar{Y}\)
- \(\hat{\sigma}_Y = \sqrt{\frac{\sum_{i=1}^n (Y_i - \overline{Y})^2}{n-1}}\)$
Sample standard error, of the sample mean, \(\bar{Y}\)
- Sample standard error, \(\hat{\sigma}_{\bar{Y}}\), is based on the sampling distribution of the sample mean \(\bar{Y}\)
- Sample standard error, \(\hat{\sigma}_{\bar{Y}}\), is the average distance between one random sample mean, \(\bar{Y}\) and the mean of the sample means, \(\bar{Y}_{\bar{Y}}\)
- note: assuming random sampling and a sufficiently large sample size, the mean of sample means, \(\bar{Y}_{\bar{Y}}\), is equal to the population mean, \(\mu_{Y}\)
- In other words, Sample standard error, \(\hat{\sigma}_{\bar{Y}}\), is the standard deviation of the sampling distribution
- Mathematically
- \(\hat{\sigma}_{\bar{Y}} = \hat{\sigma}_{Y}/\sqrt{n}\)
- equals sample standard deviation divided by the square root of sample size
Why is standard error so important?
- Example: percent of people who will vote for Obama in 2012 election (I wrote these bullets a long time ago…)
- Standard error tells us how much statistics derived from a sample are likely to diverge from population parameters
- Sample mean, \(\bar{Y}\), is best estimate of the population mean, \(\mu_Y\),, the percent of people in the population who will vote for Obama
- Standard error provides an indication of how far away each sample mean is likely to be from the population mean
- Supporse: – Standard error=10%: On average, the sample mean from each poll is likely to be 10% away from population mean – Standard error=2%: On average, the sample mean from each poll is likely to be 2% away from population mean
Do we want standard error to be large or small? Why?
- We want standard error to be small, because small standard error means that our estimates – based on sample data – are likely to be closer to the value of the population parameter of interest
- so, small standard error means more “precise” estimates
How to make standard error smaller?
- \(\hat{\sigma}_{\bar{Y}} = \hat{\sigma}_{Y}/\sqrt{n}\)
- Standard error decreases as size of your sample increases
- If you have a large sample size, the sample mean from one random sample is likely to be close to population mean – And when sample means are close to close to population mean, then standard error is small
- Example: mean income in US, with sample size of 10 versus sample size of 2,000
- Show in interactive demonstration
5 Hypothesis testing about a population mean
5.1 What and Why hypothesis testing
Quantitative research in social sciences often proceeds as follows:
- Develop a research question (which guides our research)
- Develop one (or more) testable hypothesis based on that research question
- Obtain data necessary to test the hypothesis
- Test the hypothesis by applying an appropriate statistical test to the data
Some examples of research questions co-authors and I have answered over the years:
- What is the relationship between state appropriations and nonresident enrollment at public universities (Jaquette & Curs, 2015)?
- What is the effect of nonresident enrollment growth on the number of resident students enrolled at public research universities (Curs & Jaquette, 2017)?
- What is the effect of participation in the Mexican American Studies program on the probability of high school graduation for students in the Tucscon Unified School District (Cabrera, Milem, Jaquette, & Marx, 2014)?
- Are high schools with a higher percentage of white students more likely to receive recruiting visits from university admissions officer than high schools with a lower percentage of white students (Salazar, Jaquette, & Han, forthcoming)?
For each of these journal articles, we answered the research question by developing a “testable hypothesis” and testing that hypothesis using some statistical test
Developing testable hypothesis is central to univariate statistical analysis (one variable), bivariate statistical analysis (two variables), and multivariate statistical analysis (3+ variables, usually a regression model)
Example hypotheses for univariate, bivariate, multivariate statistical analyses
- Univariate statistics (hypothesis tests about a single population mean)
- Hypothesis: the average annual cost of attendance for graduate school (tuition + fees + living expenses) is $50,000
- Bivariate statistics [hypothesis tests about comparing two population means]
- Hypothesis: the average annual cost of attendance for graduate school (tuition + fees + living expenses) at private universities is higher than public universities
- Hypothesis: the average annual cost of attendance for graduate school (tuition + fees + living expenses) at universities in urban areas is higher than universities in suburban areas
- Multivariate statistics (usually a regression model with one dependent variable, one independent variable of interest, and one or more “control” variables)
- where dependent variable (Y) = cost of attendance; independent variable of interest (X) = private or public university; control variable = level of urbanization
- Hypothesis: cost of attendance for graduate school is higher at private universities than public universities, even after controlling for level of urbanization
Why learn how to do hypothesis testing about a single population mean when this class is supposed to be about regression (and hypothesis tests about regression models)?
- You must learn the general principles/concepts about using point-estimates from sample data to test hypotheses about population parameters
- The simplest practical application of these general principles/concepts is testing hypotheses about the value of a single population mean
- The concepts/steps for hypotheses tests about a single population mean are exactly the same as those for testing hypotheses about regression models
5.2 Overview of steps in hypothesis testing
These are the general steps in hypothesis testing:
- Hypothesis
- formally state your “null” and “alternative” hypothesis
- Assumptions
- state assumptions that are relied upon by the statistical test you are using to test your hypothesis
- Test statistic
- Using some appropriate statistical analysis, calculate the “test statistic” necessary to test your hypothesis
- p-value (means probability value)
- calculate the probability of observing a test statistic as large or larger as the one you calculated
- Alpha level/rejection region and conclusion
- decide on the “alpha level,” the p-value associated with rejection of the null hypothesis
- compare the p-value you you observed to the alpha level and make a conclusion about your hypothesis test
In real research projects, do researchers always follow these exact steps? In this exact order?
- Yes, they follow these steps
- But researchers do not necessarily follow steps in this exact order
- e.g., usually, you would decide on an “alpha level” (rejection region) prior to conducting the statistical analysis
- Often, researchers will not write out each step as formally as we will ask you all to do.
- We ask you to write out each step to give you practice. Later in the quarter, you won’t have to write out each step
Example we will use to introduce steps in hypothesis testing
- The population mean cost of attendance (COA) for full-time (resident) graduate students, \(\mu_Y\), is $28,000
How we will teach you the steps in hypothesis testing in this lecture
- First, Introduce individual steps in detail, so that you develop a deep, conceptual understanding of each step
- But when thinking about an individual step in detail, it can be hard to remember its relationship to other steps and to hypothesis testing as a whole
- Second, we will do another pratical example, where we work through all steps more quickly
- so you can get a better sense of the hypothesis testing process as a whole and the relationships between steps
5.3 Hypotheses
This section presents a more formal introduction to hypotheses, focusing on univariate statistical analyses rather than bivariate or multivariate
Recall that the goal of inferential statistics is to make statements about a population of interest based on data from a representative sample from the population.
- We make a hypotheses about a population parameter (e.g., population mean of variable \(Y\) denotes \(\mu_Y\))
- Knowing the true value of the population parameter would require having data on all observations in the population
- Usually, we do not have data on the entire population
- We use sample data to test hypotheses about the population
Definition
- In statistics, a hypothesis is a declarative statement about a population.
In univariate statistical analyses, we make a hypothesis about one population paramaeter (e.g., population mean \(\mu_Y\)) from one population of interest (e.g., all “research” universities and “master’s” universities, as defined by the Carnegie Classification)
df_ipeds_pop %>% head()#> # A tibble: 6 x 38
#> instnm unitid opeid6 opeid control c15basic stabbr city zip locale region tuit_grad_res fee_grad_res tuit_grad_nres fee_grad_nres tuit_md_res fee_md_res tuit_md_nres fee_md_nres tuit_law_res fee_law_res tuit_law_nres fee_law_nres books_supplies roomboard_off oth_expense_off tuitfee_grad_res tuitfee_grad_nres tuitfee_md_res tuitfee_md_nres tuitfee_law_res tuitfee_law_nres coa_grad_res coa_grad_nres coa_md_res coa_md_nres coa_law_res coa_law_nres
#> <chr> <dbl> <chr> <chr> <dbl+lbl> <dbl+lbl> <chr+lbl> <chr> <chr> <dbl+lbl> <dbl+lbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
#> 1 Alabama A & M University 100654 001002 00100200 1 [Public] 18 [Master^s Colleges & Universities: Larger Programs] AL [Alabama] Normal 35762 12 [City: Midsize] 5 [Southeast AL AR FL GA KY LA MS NC SC TN VA WV] 10128 1414 20160 1414 NA NA NA NA NA NA NA NA 1600 9240 3090 11542 21574 NA NA NA NA 25472 35504 NA NA NA NA
#> 2 University of Alabama at Birmingham 100663 001052 00105200 1 [Public] 15 [Doctoral Universities: Highest Research Activity] AL [Alabama] Birmingham 35294-0110 12 [City: Midsize] 5 [Southeast AL AR FL GA KY LA MS NC SC TN VA WV] 8100 0 19188 0 28978 0 62714 0 NA NA NA NA 1200 12307 5555 8100 19188 28978 62714 NA NA 27162 38250 48040 81776 NA NA
#> 3 Amridge University 100690 025034 02503400 2 [Private not-for-profit] 20 [Master^s Colleges & Universities: Small Programs] AL [Alabama] Montgomery 36117-3553 12 [City: Midsize] 5 [Southeast AL AR FL GA KY LA MS NC SC TN VA WV] 11700 1300 11700 1300 NA NA NA NA NA NA NA NA 900 9600 1600 13000 13000 NA NA NA NA 25100 25100 NA NA NA NA
#> 4 University of Alabama in Huntsville 100706 001055 00105500 1 [Public] 16 [Doctoral Universities: Higher Research Activity] AL [Alabama] Huntsville 35899 12 [City: Midsize] 5 [Southeast AL AR FL GA KY LA MS NC SC TN VA WV] 10632 826 24430 826 NA NA NA NA NA NA NA NA 2120 10400 3994 11458 25256 NA NA NA NA 27972 41770 NA NA NA NA
#> 5 Alabama State University 100724 001005 00100500 1 [Public] 19 [Master^s Colleges & Universities: Medium Programs] AL [Alabama] Montgomery 36104-0271 12 [City: Midsize] 5 [Southeast AL AR FL GA KY LA MS NC SC TN VA WV] 7416 2740 14832 2740 NA NA NA NA NA NA NA NA 1600 7320 4228 10156 17572 NA NA NA NA 23304 30720 NA NA NA NA
#> 6 The University of Alabama 100751 001051 00105100 1 [Public] 16 [Doctoral Universities: Higher Research Activity] AL [Alabama] Tuscaloosa 35487-0100 12 [City: Midsize] 5 [Southeast AL AR FL GA KY LA MS NC SC TN VA WV] 10780 0 30250 0 28978 0 62714 0 23610 0 43060 0 1000 13636 4600 10780 30250 28978 62714 23610 43060 30016 49486 48214 81950 42846 62296mean(df_ipeds_pop$coa_grad_res, na.rm = TRUE)#> [1] 30318.745.3.1 null and alternative hypothesis
When developing a hypothesis for quantitative research, we always specify a null hypothesis (\(H_0\)) AND an alternative hypothesis (\(H_a\))
Null hypothesis (\(H_0\))
- In univariate statistics, a null hypothesis (\(H_0\)) is a declarative statement that the population parameter has a specific value
- (in words) \(H_0:\) the population mean cost of attendance for for full-time (resident) graduate students, \(\mu_Y\), is $28,000
- (using symbols) \(H_0: \mu_Y = \mu_{Y0} = \$28,000\)
- where \(\mu_{Y0}\) refers to the parameter value associated with the null hypothesis
- when testing a hypothesis about a single population mean, we can refer to \(\mu_{Y0}\) as the “null population mean”
Alternative hypothesis (\(H_a\))
- An alternative hypothesis (\(H_a\)) is a declarative statement that the population parameter falls in some alternative range of values as compared to the value declared by the null hypothesis
- There are two kinds of alternative hypotheses: two-sided; and one-sided
- for a given null hypothesis (\(H_0\)), there will always be one two-sided alternative hypothesis and two different one-sided hypotheses
Two-sided alternative hypothesis
- (in words) \(H_a:\) the population mean mean cost of attendance for for full-time (resident) graduate students, \(\mu_Y\), is not equal to $28,000
- (using symbols) \(H_a: \mu_Y \ne \$28,000\)
One-sided alternative hypothesis (mean is greater than \(\$28,000\))
- (in words) \(H_a:\) the population mean mean cost of attendance for for full-time (resident) graduate students, \(\mu_Y\), is greater than $28,000
- (using symbols) \(H_a: \mu_Y \gt \$28,000\)
One-sided alternative hypothesis (mean is less than \(\$28,000\))
- (in words) \(H_a:\) the population mean mean cost of attendance for for full-time (resident) graduate students, \(\mu_Y\), is less than $28,000
- (using symbols) \(H_a: \mu_Y \lt \$28,000\)
ISSUE - EQUATIONS NOT FORMATTING CORRECTLY, EVEN WHEN I CHANGE INDENTATION ASSOCIATED WITH BULLETS
Example of null and alternative hypotheses for bivariate statistical analysis
Research question:
- Is the population mean annual cost of attendance for graduate school at public universities (\(\mu_{Y_{{pub}}}\)) different from the population mean annual cost of attendance for graduate school at private universities (\(\mu_{Y_{{priv}}}\))?
Null and alternative hypotheses
- null hypothesis (\(H_0\))
- (in words): \(H_0:\) the population mean annual cost of attendance for graduate school at public universities (\(\mu_{Y_{{pub}}}\)) is the same as the population mean annual cost of attendance for graduate school at private universities (\(\mu_{Y_{{priv}}}\))
- (symbols): \(H_0: \mu_{Y_{{pub}}} = \mu_{Y_{{priv}}}\)
- Two-sided alternative hypothesis
- (in words): \(H_a:\) the population mean annual cost of attendance for graduate school at public universities (\(\mu_{Y_{{pub}}}\)) is different than the population mean annual cost of attendance for graduate school at private universities (\(\mu_{Y_{{priv}}}\))
- (symbols): \(H_a: \mu_{Y_{{pub}}} \ne \mu_{Y_{{priv}}}\)
- One-sided alternative hypothesis (\(pub < priv\))
- (in words): \(H_a:\) the population mean annual cost of attendance for graduate school at public universities (\(\mu_{Y_{{pub}}}\)) is less than than the population mean annual cost of attendance for graduate school at private universities (\(\mu_{Y_{{priv}}}\))
- (symbols): \(H_a: \mu_{Y_{{pub}}} \lt \mu_{Y_{{priv}}}\)
- note: this is the same as a one-sided hypothesis where we hypothesize \(priv > pub\)
- One-sided alternative hypothesis ($pub > priv $)
- (in words): \(H_a:\) the population mean annual cost of attendance for graduate school at public universities (\(\mu_{Y_{{pub}}}\)) is greater than than the population mean annual cost of attendance for graduate school at private universities (\(\mu_{Y_{{priv}}}\))
- (symbols): \(H_a: \mu_{Y_{{pub}}} \gt \mu_{Y_{{priv}}}\)
- note: this is the same as a one-sided hypothesis where we hypothesize \(priv < pub\)
5.3.2 Two-sided or one-sided alternative hypotheses?
In real research projects, we are not usually testing a hypothesis about a single population mean (univariate analysis). Rather, we are usually comparing population means of two different groups (bivariate analysis) or we are examining the relationship between an independent variable and the dependent variable after controlling for other variables (multivariate regression analysis)
Prior to conducting analyses, we usually have an expectation/suspicion about the result
- For most bivariate analyses, we usually suspect that one particular group is has a higher mean value than the other
- e.g., we suspect that mean cost of attendance at private universities
- this suggests a one-sided alternative hypothesis \(H_a\)
- For most multivariate analyses, we usually suspect the direction of the relationship between \(X\) and \(Y\)
- e.g., we expect that “hours spent studying” (\(X\)) has a positive relationship with “grade point average” (\(Y\)) rather than thinking "the relationship between hours spent studying (\(X\)) and grad point average (\(Y\)) does not equal zero
- this suggests a one-sided alternative hypothesis \(H_a\)
Should we specify two-sided or one-sided alternative hypothesis, \(H_a\)?
- For univariate and bivariate statistical analyses, researchers specify a two-sided alternative hypothesis more often than a one-sided alternative alternative hypothesis
- often, researchers specify a two-sided alternative hypothesis even when they strongly believe one particular group has a larger population mean then the other
- For multivariate regression analyses, researchers always specify and test two-sided alternative hypotheses
- even when they strongly believe the relationship between \(X\) and \(Y\) is positive; and even when they strongly believe the relationship between \(X\) and \(Y\) is positive
- Why this preference for two-sided alternative hypotheses in real research projects?
- two-sided alternative hypotheses are more “conservative” than one-sided alternative hypotheses;
- that is, if you specify a two-sided alternative hypothesis, \(H_a\), and reject the null hypothesis, \(H_0\), then it is necessarily true that we would have rejected then null hypothesis, \(H_0\), had we specified a one-sided alternative hypothesis, \(H_a\)
Because EDUC152 is a course about regression rather than univariate/bivariate statistics, we will always specify two-sided alternative hypotheses and ignore one-sided alternative hypotheses from here on out
- this will make learning hypothesis testing easier
5.4 Test statistic
restate null, \(H_0\), and alternative (two-sided), \(H_a\), hypothesis for our practical example
- \(H_0\)
- (in words) \(H_0:\) the population mean cost of attendance for for full-time (resident) graduate students, \(\mu_Y\), is $28,000
- (symbols) \(H_0: \mu_Y = \mu_{Y0} = \$28,000\)
- \(H_a\)
- (in words) \(H_a:\) the population mean mean cost of attendance for for full-time (resident) graduate students, \(\mu_Y\), is not equal to $28,000
- (symbols) \(H_a: \mu_Y \ne \$28,000\)
We must conduct a formal statistical test to decide whether we should reject the null hypothesis
- this is true for testing hypotheses about a population mean from a single population; testing hypotheses about whether two population means are equal; testing hypotheses about a regression coefficient, etc.
- key to understanding hypothesis testing: We conduct our test under the assumption that the null hypothesis, \(H_0\), is true
Logic of the test statistic
- What the test statistic calculates
- if the null hypothesis is true, how unlikely would it be to randomly draw the sample estimate (e.g., sample mean \(\bar{Y}\)) at least as far away from the null hypothesis value as the one we observed in our single random sample
- e.g., "if the null hypothesis is true, there is a 1.5% chance of observing a sample mean at least as far away from the null hypothesis value (\(\mu_{Y0} = \$28,000\)) as the one we observed in our single random sample
- Logic of the test statistic
- if – under the assumption that the null hypothesis is true – it would be very unlikely to observe the sample estimate we observed, then it is unlikely that the null hypothesis is true
General formula for test statistic (for pretty much any kind of hypothesis test):
\[ test\_statistic = \frac{sample\_estimate - value\_associated\_with\_H_0}{sample\_standard\_error}\] Formula for test statistic about a single population mean
- in words:
- test-statistic \(t\) equals difference between the sample mean \(\bar{Y}\) and the population mean associated with the null hypothesis \(\mu_{Y0}\) divided by the sample standard error of the sample mean \(\hat{\sigma}_{\bar{Y}}\)
- equation:
- \(t = \frac{\bar{Y} - \mu_{Y0}}{\hat{\sigma}_{\bar{Y}}}\)
- where:
- \(\hat{\sigma}_{Y}\) refers to sample standard deviation of variable \(Y\)
- \(n\) refers to sample size
- sample standard error of the sample mean \(= \hat{\sigma}_{\bar{Y}} = \frac{\hat{\sigma}_{Y}}{\sqrt{n}}\)
Calculating t-test statistic for our practical example
\(H_0: \mu_Y = \mu_{Y0} = \$28,000\) ; \(H_a: \mu_Y \ne \$28,000\)
Calculate components of t-test (using functions and by hand)
# sample size
length(df_ipeds_sample$coa_grad_res) # assuming no missing observations#> [1] 200 df_ipeds_sample %>% summarize(n_non_miss = sum(!(is.na(coa_grad_res)))) # count only number of non-missing#> # A tibble: 1 x 1
#> n_non_miss
#> <int>
#> 1 200# sample mean of coa_grad_res
mean(df_ipeds_sample$coa_grad_res, na.rm = TRUE) # using function#> [1] 30002.74# sample standard deviation of coa_grad_res
sd(df_ipeds_sample$coa_grad_res, na.rm = TRUE)#> [1] 11457.01# sample standard error of sample mean of coa_grad_res = std_dev/sqrt(n)
sd(df_ipeds_sample$coa_grad_res, na.rm = TRUE)/sqrt(length(df_ipeds_sample$coa_grad_res))#> [1] 810.1332Components of t-test:
- sample size, \(n\) = 200
- sample mean, \(\bar{Y}\) = 30002.745
- Population mean associated with \(H_0\), \(\mu_{Y0} = \$28,000\)
- sample standard deviation, \(\hat{\sigma}_{Y}\) = 11457.014
- sample standard error of the sample mean, \(\hat{\sigma}_{\bar{Y}}\) = 810.1332
Calculating t-test
\[t = \frac{\bar{Y} - \mu_{Y0}}{\hat{\sigma}_{\bar{Y}}} = \frac{30002.74 - 28000}{810.1332} = 2.4721\]
t.test() function
- what
t.test()does- “Performs one and two sample t-tests on vectors of data”
- we use one sample t-test to test hypothesis about single population mean
- later, we will use two-samle t-test to test hypotheses about whether two population means are equal
- syntax:
t.test(x, y = NULL, alternative = c("two.sided", "less", "greater"), mu = 0, paired = FALSE, var.equal = FALSE, conf.level = 0.95, ...)
- selected arguments
x: vector (variable) you want to calculate t-test foralternative: whether you want two-sided or one-sided alternative hypothesis (default istwo.sided)mu: value associated with null hypothesis (default is0)
Calculating t-test value (using function and by hand)
# t-statistic = (sample_mean - mu_H_0)/(sample std err)
# using function
#?t.test # to see help file for function
t.test(x = df_ipeds_sample$coa_grad_res, mu = 28000)#>
#> One Sample t-test
#>
#> data: df_ipeds_sample$coa_grad_res
#> t = 2.4721, df = 199, p-value = 0.01427
#> alternative hypothesis: true mean is not equal to 28000
#> 95 percent confidence interval:
#> 28405.20 31600.29
#> sample estimates:
#> mean of x
#> 30002.74# by hand
(mean(df_ipeds_sample$coa_grad_res, na.rm = TRUE) - 28000)/(sd(df_ipeds_sample$coa_grad_res, na.rm = TRUE)/sqrt(length(df_ipeds_sample$coa_grad_res)))#> [1] 2.4721185.4.1 conceptual understanding of test statistic (MOST IMPORTANT!)
The test statistic refers to a sampling distribution not the distribution of your single sample
- in particular, the test statistic refers to the sampling distribution under the assumption that the null hypothesis is true, \(H_0: \mu_Y = \mu_{Y0}\)
Recall core ideas of sampling distribution (we will refer to sampling distribution of sample mean, \(\bar{Y}\))
- The sampling distribution of the sampling mean is a distribution where each observation is a sample mean, \(\bar{Y}\), from one random sample taken from the population
- The sampling distribution shows how the value of the sample mean varies from sample to sample
- Standard error (i.e., the standard deviation of the sampling distribution), \(\hat{\sigma}_{\bar{Y}}\) is the average distance between one random sample mean, \(\bar{Y}\), and the mean of sample means \(\bar{Y}_{\bar{Y}}\)
- (Note: we refer to sample standard error of the sample mean \(\hat{\sigma}_{\bar{Y}} = \frac{\hat{\sigma}_{Y}}{\sqrt{n}}\), which can be calculated from sample data, rather than population standard error of the sample mean \(\sigma_{\bar{Y}}\), which is a population parameter that is only known if we have data on the entire population)
- Drawing from the central limit theorem, we know that sampling distributions will always be normally distributed so long as sample size is not small
- Therefore, sampling distributions follow the empirical rule:
- 68% of observations within one standard error
- 95% of observations wtithin two standard erros
- 99% of observations within three standard errors
Here we visually stack the following for the variable coa_grad_res:
- the population distribution; distribution from a single random sample; and sampling distribution of sample mean
plot_distribution(df_ipeds_pop, 'coa_grad_res', plot_title = 'Population distribution') +
plot_distribution(df_ipeds_sample, 'coa_grad_res', plot_title = 'Single sample distribution') +
plot_distribution(df_ipeds_pop, 'coa_grad_res', sampling_dist = T,
plot_title = 'Sampling distribution') +
plot_layout(ncol = 1)
Usually we cannot know the sampling distribution because we do not have data on the entire population; we only have data on our single random sample
However, hypothesis testing is not based on the true sampling distribution of the sample mean. It is based on the sampling distribution under the assumption that the null hypothesis is correct
- Thanks to the central limit theorem, we have a pretty good idea of the sampling distribution assuming \(H_0\) is true, even when we only have a single random sample!
Here we visually stack the following for the variable coa_grad_res:
- the population distribution; distribution from a single random sample; and sampling distribution of assuming that \(H_0\) is true
plot_distribution(df_ipeds_pop, 'coa_grad_res', plot_title = 'Population distribution') +
plot_distribution(df_ipeds_sample, 'coa_grad_res', plot_title = 'Single sample distribution') +
plot_t_distribution(df_ipeds_sample, 'coa_grad_res', mu = 28000,shade_rejection = F, shade_pval = T,
plot_title = 'Sampling distribution, assuming H_0', stacked = T) +
plot_layout(ncol = 1)
The t-test statistic is the distance between the hypothesized \(H_0\) value and the observed sample estimate value \(\bar{Y}\) scaled in terms of standard errors
- e.g., it would be unlikely to observe a t-value of greater than
2or a t-value less than-2because we know (from empirical rule and central limit theorem) that 95% of observations fall within two standard deviations of the mean for a normally distributed variable
5.5 p-value
“p-value” refers to the probability-value associated with the t-value from your test statistic
Definition
- Under the assumption that \(H_0\) is true, the p-value is the probability of observing a sample estimate (and its associated test-statistic) that is at least as far away from the null hypothesis value \(\mu_{Y0}\) as the one we observed
A small p-value means that it would be unusual to find the sample estimate we observed if the null hypothesis \(H_0\) is find the observed data if 𝐻0were true.
Calculating p-value For a two-sided alternative hypothesis (\(H_a: \mu_Y \ne \mu_{Y0}\))
- let \(t\) be the value of your t-test
- let \(p\) by the p-value associated with \(t\)
- let \(Pr(obs>t)\) is the probability of an observation having a higher value of \(t\) than the one you observed
- \(p = Pr(obs > t) + Pr(obs< -t)\)
- Because the sampling distribution is symmetric (because it is normally distributed):
- \(Pr(obs >t) = Pr(obs < - t)\)
- therefore, for a two-sided alternative hypothesis
- \(p = 2*Pr(obs > t)\)
Let’s calculate and visualize p-value for a couple different hypothesized values of the population mean \(\mu_{Y0}\) for the variable coa_grad_res (full-time, resident grad school cost of attendance) from the data frame df_ipeds_sample
\(H_0: \mu_Y = \mu_{Y0} = \$29,000\) and \(H_a: \mu_Y \ne \$29,000\)
- Sample mean, ${Y} = $ 30002.745
- \(t =\) 1.24
- p-value \(= Pr(obs > t) + Pr(obs< -t) =\) 0.217
- \(Pr(obs>t)=\) 0.1086
- \(Pr(obs<-t)=\) 0.1086
- below code chunk runs t-test and plots t-value against sampling distribution assuming \(H_0\) is true
mean(x = df_ipeds_sample$coa_grad_res)#> [1] 30002.74t.test(x = df_ipeds_sample$coa_grad_res, mu = 29000)#>
#> One Sample t-test
#>
#> data: df_ipeds_sample$coa_grad_res
#> t = 1.2378, df = 199, p-value = 0.2173
#> alternative hypothesis: true mean is not equal to 29000
#> 95 percent confidence interval:
#> 28405.20 31600.29
#> sample estimates:
#> mean of x
#> 30002.74plot_t_distribution(df_ipeds_sample, 'coa_grad_res', mu = 29000,shade_rejection = F, shade_pval = T)
\(H_0: \mu_Y = \mu_{Y0} = \$28,000\) and \(H_a: \mu_Y \ne \$28,000\)
- Sample mean, ${Y} = $ 30002.745
- \(t =\) 2.47
- p-value \(= Pr(obs > t) + Pr(obs< -t) =\) 0.014
- \(Pr(obs>t)=\) 0.0071
- \(Pr(obs<-t)=\) 0.0071
- below code chunk runs t-test and plots t-value against sampling distribution assuming \(H_0\) is true
mean(x = df_ipeds_sample$coa_grad_res)#> [1] 30002.74t.test(x = df_ipeds_sample$coa_grad_res, mu = 28000)#>
#> One Sample t-test
#>
#> data: df_ipeds_sample$coa_grad_res
#> t = 2.4721, df = 199, p-value = 0.01427
#> alternative hypothesis: true mean is not equal to 28000
#> 95 percent confidence interval:
#> 28405.20 31600.29
#> sample estimates:
#> mean of x
#> 30002.74plot_t_distribution(df_ipeds_sample, 'coa_grad_res', mu = 28000,shade_rejection = F, shade_pval = T)
\(H_0: \mu_Y = \mu_{Y0} = \$31,500\) and \(H_a: \mu_Y \ne \$31,500\)
- Sample mean, \(\bar{Y} =\) 30002.745
- \(t =\) -1.85
- p-value \(= Pr(obs > t) + Pr(obs< -t) =\) 0.066
- \(Pr(obs>t)=\) 0.033
- \(Pr(obs<-t)=\) 0.033
- below code chunk runs t-test and plots t-value against sampling distribution assuming \(H_0\) is true
mean(x = df_ipeds_sample$coa_grad_res)#> [1] 30002.74t.test(x = df_ipeds_sample$coa_grad_res, mu = 31500)#>
#> One Sample t-test
#>
#> data: df_ipeds_sample$coa_grad_res
#> t = -1.8482, df = 199, p-value = 0.06606
#> alternative hypothesis: true mean is not equal to 31500
#> 95 percent confidence interval:
#> 28405.20 31600.29
#> sample estimates:
#> mean of x
#> 30002.74plot_t_distribution(df_ipeds_sample, 'coa_grad_res', mu = 31500)
5.6 Rejection region and conclusion
5.6.1 Alpha level (rejection region)
\(\alpha\) level (referred to as “alpha level” or “rejection region”)
- Definition: \(\alpha\) level is a number such that we reject the null hypothesis \(H_0\) the observed p-value is less than or equal to the alpha level.
In practice, the most common alpha level \(\alpha\) is .05
- e.g., if we choose \(\alpha= .05\) and our t-statistic is associated with a p-value of .02, then we reject \(H_0\); if we choose \(\alpha= .05\) and our t-statistic is associated with a p-value of .07, then we do not reject \(H_0\)
- sometimes researchers choose an alpha level of
.10but usually this is viewed as not sufficiently strong threshold to reject \(H_0\)
Usually, you define alpha level prior to running analyses
- when researchers define alpha level (rejection region) after running analyses, there may be a temptation to choose an alpha level that allows them to reject \(H_0\)
To show how alpha level is used in practice, we’ll test the null hypothesis that population mean grad resident cost of attendance is 28,500, initially using an alpha level of .05
- \(H_0: \mu_Y = \mu_{Y0} = \$28,500\) and \(H_a: \mu_Y \ne \$28,500\)
- Sample mean, \(\bar{Y} =\) 30002.745
- \(t =\) 1.85
- p-value \(= Pr(obs > t) + Pr(obs< -t) =\) 0.065
- \(Pr(obs>t)=\) 0.0325
- \(Pr(obs<-t)=\) 0.0325
Note that the t.test() function doesn’t have an argument that let’s you specify the alpha level (rejection region); rather, the idea is that you choose the alpha level and then compare that to the p-value calculated by t.test()
t.test(x = df_ipeds_sample$coa_grad_res, mu = 28500)#>
#> One Sample t-test
#>
#> data: df_ipeds_sample$coa_grad_res
#> t = 1.8549, df = 199, p-value = 0.06508
#> alternative hypothesis: true mean is not equal to 28500
#> 95 percent confidence interval:
#> 28405.20 31600.29
#> sample estimates:
#> mean of x
#> 30002.74Note that the user-defined plot_t_distribution() includes an optional argument that allows you to specify the alpha-level
- syntax (including default argument values for arguments with defaults)
plot_t_distributionfunction(data_df, data_var, group_var = NULL, group_cat = NULL, mu = 0, alpha = 0.05, alternative = 'two.sided', plot_title = '', shade_rejection = T, shade_pval = T, stacked = F)- So we can set the alpha level with the argument
alpha, which has the default value of0.05
Below, we run plot_t_distribution(), manually setting the alpha argument to the the default value of 0.05
- the blue dotted line denotes the t-value from our test statistic
- \(t =\) 1.85
- the red dotted lines denote the t-value associated with the chosen alpha level (rejection region)
- the t-value associated with the alpha level is referred to as the “critical value”
- the t-statistic of 1.85 is less than the critical value of 1.97, indicating that we won’t reject \(H_0\)
- the shaded region denotes the probability associated with the chosen alpha level
- given our choice of
alpha = .05, the shaded region represents .05, that is the percent of all observations that lie in the rejection region - our p-value of 0.217 is greater than our alpha level of .05, indicating that we won’t reject \(H_0\)
- given our choice of
- Since our p-value is greater than our alpha level, we do not reject \(H_0\)
plot_t_distribution(df_ipeds_sample, 'coa_grad_res', mu = 28500, alpha = .05)
How would our conclusion change if we chose an alpha level of .10?
- \(H_0: \mu_Y = \mu_{Y0} = \$28,500\) and \(H_a: \mu_Y \ne \$28,500\)
- Sample mean, \(\bar{Y} =\) 30002.745
- \(t =\) 1.85
- p-value \(= Pr(obs > t) + Pr(obs< -t) =\) 0.065
- \(Pr(obs>t)=\) 0.0325
- \(Pr(obs<-t)=\) 0.0325
- Since our p-value of 0.065 is less than the alpha level of .10, we reject \(H_0\)
plot_t_distribution(df_ipeds_sample, 'coa_grad_res', mu = 28500, alpha = .10)
5.6.2 Conclusion
Last step is to make a conclusion about your hypothesis based on comparing the p-value you observed to the alpha level
How to state the conclusion
- if p-value is greater than the alpha level \(\alpha\) (e.g., \(p-val =.06, \alpha = .05\))
- conclusion: do not reject \(H_0\)
- longer version: “We do not have sufficient evidence to reject the null hypothesis, \(H_0\), that population mean \(\mu_Y\) is equal to \(\mu_{Y0}\)”
- do not say "we accept \(H_0\) because we don’t know the exact value of the population parameter
- if p-value is less than the alpha level \(\alpha\) (e.g., \(p-val =.04, \alpha = .05\))
- conclusion: we reject \(H_0\)
- can also say “we accept \(H_a\), but people usually say”reject \(H_0\)
- longer version: “we reject the null hypothesis \(H_0\), that population mean \(\mu_Y\) is equal to \(\mu_{Y0}\)”
- conclusion: we reject \(H_0\)
This example:
- \(H_0: \mu_Y = \mu_{Y0} = \$28,500\) and \(H_a: \mu_Y \ne \$28,500\) and alpha level \(\alpha = .05\)
- p-value \(= Pr(obs > t) + Pr(obs< -t) =\) 0.065
- conclusion:
- do not reject \(H_0\)
- We do not have sufficient evidence to reject the null hypothesis, \(H_0\), that population mean full-time resident graduate cost of attendance is equal to
28,500
plot_t_distribution(df_ipeds_sample, 'coa_grad_res', mu = 28500, alpha = .05,
shade_rejection = TRUE, shade_pval = FALSE)
5.7 Assumptions
All statistical tests (based on some statistical analysis) depends on “assumptions”
- if the researcher is confident that the assumptions have been satisfied, then the researcher can make inferences about the population parameter by applying the relevant statistical analysis/test to sample data
- if we are concerned that one or more assumptions have not been satisfied, then the researcher should not make inferences about the population parameter
Assumptions necessary for testing hypotheses about a population mean
- sample is a random sample from population
- population distribution of variable is normal
“Robust”
- A statistical method is robust with respect to a particular assumption, when it performs adequately even when that assumption is violated
Hyptohesis tests about a population mean is “robust” to the normal distribution assumption
- Statisticians have shown that hypothesis tests about population means are robust against violations of normal population assumption, especially when sample size > 30
- Why is hypothesis test about population means robust to normal population assumption? Because of central limit theorem
Central limit theorem:
- when sample size is large, the sampling distribution of the sample mean, 𝑦 , is approximately normal, even if the population distribution of the variable is not normal
- If population distribution is normal then sampling distribution is normal for any sample size
- If population distribution is not normal, sample size of about 30 is sufficient
Hypothesis test about population means is not robust to violations of random sampling
- i.e., if you take a non-random sample from the population, you cannot make good predictions about the population
5.8 Hypothesis test example, all steps
Research question:
- What is the population mean price of full-time nonresident graduate tuition + fees? [variable =
tuitfee_grad_nres]
Let’s imagine we want to test whether the population mean, \(\mu_Y = \$17,000\), using a two-sided alternative hypothesis and an alpha level of .05
State null and alternative hypotheses
- Null hypothesis, \(H_0\)
- \(H_0: \mu_Y = \mu_{Y0} = \$17,000\)
- \(H_0:\) population mean price of full-time nonresident graduate tuition + fees is $17,000
- Alternative hypothesis, \(H_a\)
- \(H_0: \mu_Y \ne \$17,000\)
Test statistic and p-value
- \(t = \frac{\bar{Y} - \mu_{Y0}}{\hat{\sigma}_{\bar{Y}}}\)
- where:
- \(\hat{\sigma}_{Y}\) refers to sample standard deviation of variable \(Y\)
- \(n\) refers to sample size
- sample standard error of the sample mean $ = _{{Y}} = $
Components of t-test
- sample size, \(n\) = 200
- sample mean, \(\bar{Y}\) = 18616.63
- Population mean associated with \(H_0\), \(\mu_{Y0} = \$17,000\)
- sample standard deviation, \(\hat{\sigma}_{Y}\) = 9976.3157
- sample standard error of the sample mean, \(\hat{\sigma}_{\bar{Y}}\) = 705.432
Calculating t-test p-value using t.test()
- \(t =\) 2.29
- p-value \(= Pr(obs > t) + Pr(obs< -t) =\) 0.023
- \(Pr(obs>t)=\) 0.0115
- \(Pr(obs<-t)=\) 0.0115
- below code chunk runs t-test and plots t-value against sampling distribution assuming \(H_0\) is true, using alpha of .05
t.test(df_ipeds_sample$tuitfee_grad_nres, mu = 17000)#>
#> One Sample t-test
#>
#> data: df_ipeds_sample$tuitfee_grad_nres
#> t = 2.2917, df = 199, p-value = 0.02297
#> alternative hypothesis: true mean is not equal to 17000
#> 95 percent confidence interval:
#> 17225.55 20007.71
#> sample estimates:
#> mean of x
#> 18616.63plot_t_distribution(df_ipeds_sample, 'tuitfee_grad_nres', mu = 17000, alpha = .05,
shade_rejection = TRUE, shade_pval = FALSE)
Conclusion
- Because the p-value of 0.023 is less than the alpha level of
.05, we reject \(H_0\) - we reject the null hypothesis \(H_0\), population mean price of full-time nonresident graduate tuition + fees, \(\mu_Y\), is equal to \(17,000\)
- We can also say that \(\mu_Y\) is greater than \(17,000\)
Finally, we usually don’t have all data on the population. But since we do for IPEDS, we can plot:
- the population distribution (usually unknown)
- on top of the distribution from our single random sample
- on top of the sampling distribution (usually unknown)
- on top of the sampling distribution assuming \(H_0\) is true
plot_distribution(df_ipeds_pop, 'tuitfee_grad_nres', plot_title = 'Population distribution') +
plot_distribution(df_ipeds_sample, 'tuitfee_grad_nres', plot_title = 'Single sample distribution') +
plot_distribution(df_ipeds_pop, 'tuitfee_grad_nres', sampling_dist = T, plot_title = 'True Sampling distribution') +
plot_t_distribution(df_ipeds_sample, 'tuitfee_grad_nres', mu = 17000, plot_title = 'Sampling distribution, assuming H_0', stacked = T) +
plot_layout(ncol = 1)
6 Fundamentals of causal inference
text
7 Hypothesis testing about two populations
text
References
Cabrera, N. L., Milem, J. F., Jaquette, O., & Marx, R. (2014). Missing the (student achievement) forest for all the (political) trees: Empiricism and the mexican american studies controversy in tucson. American Educational Research Journal, 51(6), 1084–1118.
Curs, B. R., & Jaquette, O. (2017). Crowded out? The effect of nonresident enrollment on resident access to public research universities. Educational Evaluation and Policy Analysis, 39(4), 644–669. https://doi.org/10.3102/0162373717704719
Jaquette, O., & Curs, B. R. (2015). Creating the out-of-state university: Do public universities increase nonresident freshman enrollment in response to declining state appropriations? Research in Higher Education, 56(6), 535–565. https://doi.org/10.1007/s11162-015-9362-2
Salazar, K., Jaquette, O., & Han, C. (forthcoming). Coming soon to a neighborhood near you? Off-campus recruiting by public research universities. American Educational Research Journal. Forthcoming.